Correlation lengths in hydrodynamic models of active nematics
\makeFNbottom\sectionfont\subsectionfont\subsubsectionfont\setstretch1.125 0.8cm \titlespacing*
0pt4pt4pt \titlespacing*
0pt15pt1pt
2 Introduction
Active systems are continuously driven out of equilibrium by energy injected at the local scale, resulting in collective motion at the large scale. Examples include bacterial colonies, in vitro extracts of cytoskeletal filaments and associated motor proteins, and monolayers of vibrated granular matter ^{1}. Much interest has recently focused on active nematics  systems composed of active units with headtail symmetry that can order into states with nematic liquid crystalline order, with the activity giving rise to a rich variety of collective phenomena. These include spontaneous laminar flow ^{2, 3, 4}, large density fluctuations ^{5, 6, 7}, pattern formation ^{8, 9}, spontaneous unbinding of topological defects ^{10, 11, 12, 13}, and active turbulence ^{14, 15, 16, 17, 18}.
Theoretical interest in active nematics has been fueled by the observation of spontaneously flowing and turbulent states in suspensions of microtubulekinesin bundles confined at an oilwater interface ^{19, 20, 21}. Other recent experimental realizations were obtained by immersing living swimming bacteria (specifically, E. coli) in a lyotropic liquid crystal ^{22} and by plating dense layers of fibroblasts on substrates ^{23}. Microtubule bundles and E. coli exert active extensile force dipoles on their surroundings, while the fibroblasts exert contractile force dipoles. In both systems active stresses couple to orientational order and induce flows and defect unbinding, with qualitatively similar nonequilibrium dynamics, albeit on very different time scales.
All these experiments have indicated that the appearance of spatiotemporal chaos in active nematics is accompanied by the proliferation of topological defects, which may mediate the onset of turbulence even at negligible values of the Reynolds number ^{11, 14}. In active fluids, distortions of the local orientational order produce local shear flows that enhance the orientational deformation, ultimately leading to the unbinding of defect pairs. In nematics, these consist of pairs of disclinations  specific distortions of the orientational order that are the signature of the underlying broken symmetry. It was additionally shown ^{11} that the dynamics of topological defects in active systems depends on the nature of the forcing, i.e., whether the force dipole is contractile or extensile. This result provides a useful criterion for inferring the nature of active forces in systems with nematic symmetry ^{11, 16, 24}.
In spite of much theoretical work, discrepancies still exist in the literature over the nature of the characteristic length scales that control the spontaneous proliferation and annihilation of topological defects, and the resulting dynamics in the socalled turbulent state. In particular, the dependence of such length scales on the strength of the active forcing, , remains unclear. A recent numerical study by Giomi^{17} examined the statistics of the activitydriven turbulent phase in twodimensional nematic films by measuring the distribution of vortex sizes for a selection of activities. This work provided evidence that the key physics is determined by a single active length scale, , proportional to . In contrast, in a closely related work, but on a different continuum model of a quasi nematic, Thampi et al. performed a detailed study that measured several orientational and hydrodynamical correlation lengths, suggesting that the length scale of structure in the fluid instead scales as .
An important aim of this work is to provide a unified understanding of these previous, apparently conflicting reports. In order to do so, we consider two different but related models of active nematic liquidcrystal hydrodynamics that are commonly used in the literature. The dynamics of these models is compared using two independently developed numerical codes. By varying the key dimensionless parameters over several decades, we obtain data to support the conjecture that, in both models, the mean defect spacing in the regime of full developed active turbulence is set by the length scale , defined by the balance of active and elastic stresses. Here is the magnitude of the active stress (commonly referred to as the activity) and parametrizes the free energy penalty that results from spatial variations in the director field ^{3, 4, 1, 13, 17, 25}. We show that this result holds for both extensile and contractile systems, in both the flow aligning and the flow tumbling regimes. This active length scale also controls the onset of spontaneous laminar flow in an active film, a phenomenon that has been referred to in the literature as the spontaneous flow instability ^{2}. Our study provides the first explicit demonstration that distinct constitutive models produce the same emergent length scale, i.e., they both produce quantitatively consistent scaling relations. We also demonstrate a regime of less highly developed turbulence in which a weaker scaling appears consistent with our numerical data.
In many experimental realizations active nematics are confined to quasi twodimensional geometries, e.g., on the surface of lipid vesicles ^{20}, in flattened waterinoil droplets ^{19}, in thinfilms ^{26}, or squeezed between parallel glass plates ^{22}. It is worth noting that the presence of confining walls (with noslip boundary conditions) in the last of these modifies hydrodynamic interactions between active particles ^{27}; of the above experiments, freestanding films may then provide the closest experimental realisation of our system. Clearly, in any numerical study, it is important to define carefully the considered dimensionality. In what follows we denote by the number of dimensions in which the relevant fields (nematic order parameter tensor, fluid velocity, etc.) are allowed to vary; and separately by the number of dimensions in which the nematic director is allowed to develop nonzero components. We shall perform two different studies. In the first we take a strictly twodimensional model of an active nematic sheet, in which the order parameter tensor is allowed to develop nonzero components only in the directions (); and physical quantities are likewise allowed to vary only in the plane (). In the second study we consider a threedimensional nematic () but in which all quantities are nonetheless still assumed to be spatially homogeneous in the direction of the layer thickness (). In the latter case the director can in principle point out of the simulated plane, and indeed this effect has been reported in a previous numerical study of active nematics in a cylindrical capillary ^{28}. In that work confinement in the plane of simulation plays an important role whereas our study deliberately focuses on the bulk dynamics where we do not observe outofplane motion. Finally it we note that accurately resolving turbulent active hydrodynamics is computationally demanding, especially when performing large sweeps of parameter space. This necessarily restricts us to a 2D study, as is the case with many other studies of active turbulence ^{3, 14, 17}.
The paper is structured as follows. In Section 3 we define the equations of motion for both models, and outline the parameter ranges that we explore for each. In Section 4, we define the observable length scales that can be used to characterize the fluid structure, and discuss the physical reasoning behind their definition. The results of our study are presented in Section 5. In Section 6 we provide a comparison with other works and offer our conclusions.
3 Models
In spatial dimensions we consider an incompressible uniaxial active nematic liquid crystal with a director that can orient in dimensions, with in our two respective studies. The nematic orientational order is parametrized by a symmetric and traceless tensor field , where is the order parameter magnitude and the director is a headless unit vector that characterizes the direction of broken orientational symmetry. The nematic is embedded in an incompressible fluid of constant density, , and constant viscosity, . The fluid velocity field is denoted by . The associated pressure field is determined by the incompressibility condition .
The equations of motion for an active nematic are derived from the wellknown hydrodynamic equations for a passive liquidcrystal ^{25}
(1)  
(2) 
where is the material derivative and is a rotational viscosity. Here and denote the symmetric and antisymmetric parts of the rate of strain tensor , respectively, with and . For other tensors the transpose, symmetric, antisymmetric and traceless parts are denoted by the superscripts , , and , respectively. For example, .
The relaxation dynamics of the alignment tensor in Eq. 2 is governed by the molecular field, , in which the Landaude Gennes free energy ^{25}, , is the sum of contributions from a bulk free energy density
(3) 
and the distortion free energy density
(4) 
Here and determine the bulk and distortion energy density scales respectively. For simplicity we have adopted the oneelastic constant approximation in the distortion free energy (Eq. 4).
The tensor is defined for an arbitrary tensor as
(5) 
where
(6) 
Here is the LeslieEricksen flow aligning parameter, which specifies how the nematic director responds to a shear flow: corresponds to flowaligning nematics and corresponds to the flowtumbling regime. (See Appendix A.)
The stress tensor in Eq. 1 is the sum of passive liquidcrystal and active contributions, . The passive part of the stress tensor is given by
(7) 
In an active nematic there is an additional active stress contribution that arises from the dipolar forces exerted by active particles on their environment. This active stress is , where describes contractile stresses and extensile stresses. In the passive limit , the equations just described reduce to those of a passive liquidcrystal.
So far, the model that we have presented encompasses both of the numerical studies performed. We now outline the specific choices for parameter values and dimensionality made in each of the two numerical studies separately, and discuss how these choices affects the form of the equations. The main difference between the two variations will be the presence of higher order coupling terms in Model II, both in the liquidcrystal stress and in the coupling between orientational order and velocity gradients.
Model I (MI): This model describes a dimensional nematic sheet with a dimensional nematic order parameter. In this case the symmetric second rank tensor has only two independent components and identically. In Model I the meanfield free energy (Eq. 3) has coefficients and , where is a dimensionless parameter that controls the continuous transition from an isotropic to a nematic state, with the transition occurring at . The second term in Eq. 5 is also identically zero; the third term is of a higher order in and can safely be neglected ^{29}, so that . We also exclude the last term in Eq. 7 in Model I. We assume a constant density for which the Reynolds number , where the velocity scale is
Model II (MII): This model considers a dimensional layer of nematic liquid crystal described by the full Landau free energy given in Eq. 3, thereby in principle allowing the director to explore all dimensions. However it still neglects all spatial variations in the direction of the layer thickness, so taking as noted above. In this case the free energy in Eq. 3 sets , and , yielding a first order isotropicnematic transition at . In the following we choose which places us at the spinodal stability limit of the isotropic phase and well within the nematic state, with . According to Eq. 6, the system will be in the flowaligning regime for and in the flowtumbling regime for . In all simulations using Model II we have fixed , corresponding to a flowaligning system. We consider only extensile systems with this model, i.e., values of . Finally, in MII we take the limit of zero Reynolds number by setting .
parameter  description  dimensions 

activity  
Frank constant  ( in MI)  
energy density scale  ( in MI, MII)  
rotational viscosity  ( in MII)  
alignment param.  
IN control param.  
solvent viscosity  ( in MI)  
solvent density  
box size  ( in MII) 
The full list of nine parameters (for both models) is given in Table 1. We are free to choose units of mass , length and time , or equivalently of stress , length and time , and we have noted in Table 1 which quantities we chose to set equal to unity in each of the two studies. This leaves six dimensionless groupings that we summarize in Table 2, three of which are fixed throughout. Therefore even though we choose our units differently in the two different simulation studies, all results are presented and compared in a consistent adimensional way between the two models. We choose parameters that produce flowtumbling behaviour in Model I and flowaligning in Model II. Note that due to differences in parameter selections, the linear instability thresholds in the two models differ by a factor ( in Model I ^{13} and in Model II ^{30}). Accordingly, the onset of the turbulent regime in each model is separated by a similar factor, requiring us to explore different ranges of dimensionless activity, as noted in Table 2. In particular, in Model II we explore the transition from small to large activities, whereas Model I focuses on larger activities still (i.e., deeper into the regime of fully developed active turbulence).
parameter  description  MI value  MII value 

varied parameters  
dimensionless activity  
ratio of micro to macroscopic length scales  
ratio of viscosities  0.567  
fixed parameters  
alignment param.  0.7  
IN control param.  2  3  
Reynolds number  0 
3.1 Numerical details
In order to demonstrate the robustness of our results with respect to numerical implementation, we use two independent codes (one for each model). In each case we perform simulations in a square box of side with biperiodic boundary conditions. The dynamics in Model I is timeintegrated numerically on a square grid of points using a fourth order RungeKutta method, with a timestep . Gradients of are computed using a finite difference scheme. In Model II is integrated numerically using a Euler timestepping scheme of timestep in the range on a grid of points (dependent on the magnitude of activity) and gradients of are treated using a semiimplicit Fourier method. In Model II, the velocities are determined instantaneously from the force balance equation, which we solve in Fourier space using a stream function formulation. In Model I, we integrate the NavierStokes equation (Eq. 1) with the same scheme used for the order parameter equation to obtain the velocity at every time step. We have verified that our results are quantitatively unchanged upon decreasing the timestep or grid spacing. Both simulations were initialised with a uniform director field orientated within the plane, subject to a small sinusoidal perturbation of magnitude .
4 Characteristic length scales
Irrespective of the specific details of the model used, we expect the resulting dynamics of the active nematic to be controlled by the interplay of key length and time scales that govern the basic physics.
4.1 Model length scales
An inspection of the hydrodynamic equations and model geometry reveals three underlying length scales. The first is simply the system size . The second arises from balancing the bulk and elasticdistortion free energy terms in Eqs. 3 and 4, to obtain the equilibrium nematic persistence length, which deep in the nematic state is given by
(8) 
This is the length scale over which spatial correlations in the nematic field decay deep in the nematic phase, where it is proportional to the defect core radius. The third lengthscale arises by balancing the elastic stress associated with a deformation over a length with the active stress scale , to give the active length scale
(9) 
To guarantee that any physics on these lengthscales is not contaminated by finite size effects, we focus on the regime in which and .
Alternatively, from a dynamical viewpoint one might consider the system to be controlled by two timescales: the passive structural relaxation time , which controls the relaxation of a distortion to the nematic order on a length scale , and the active time scale , which controls the relative rates of injection of active stresses and stress decay via viscous dissipation. The length scale that results when these timescales are equated is then
(10) 
4.2 Emergent length scales
The length scales discussed above were motivated by simple dimensional analysis of the model parameters and flow geometry. In our numerical simulations, we find that (for a high enough level of activity) an initially homogeneous state gives way to a spatiotemporally complicated state with defects in the nematic director field, and associated local flows in the velocity field, as found earlier by several authors ^{11, 14} and shown in the snapshots of Figs. 2 and 3. An important aim of the present work is to elucidate how the length scales associated with these emergent structures depend on the underlying model length scales just discussed. We denote these emergent length scales by the common symbol , but in fact there are multiple possible scales that we might choose to characterize the spatiotemporal dynamics, as we now describe.

Mean defect separation : We define the mean defect separation
(11) where is the areal density of defects, calculated by adapting the defect tracking method of Ref. ^{31}.

Director correlation length : The normalised director correlation function defined as
(12) This characterizes the probability that two director orientations a distance apart are the same (respecting the fact that are equivalent for a nematic). Here and throughout, the angular brackets indicate an average over space and time. We then choose to be the length at which , as in Fig. 1 (inset).

Velocity correlation length : Analogously, the velocity correlation function,
(13) defines the velocity correlation length according to .

Vorticity correlation length : Finally, we define the correlation function for the local vorticity, , as
(14) and define the vorticity correlation length by .
4.3 Scaling hypothesis
Simple dimensional analysis based on the model length scales discussed in Sec. 4.1 suggests that the length scales of Sec. 4.2 characterizing the emergent structures in the fluid (whether or ) should obey a simple scaling relation of the form
(15) 
where is a general scaling function.
Previous simulation studies^{3, 32, 30, 17, 15} have shown that all characteristic length scales, denoted generically by , decrease with increasing activity . At low activity, typically just a few defects are seen in the simulation box, as in the snapshots in Figs. 2c and 3b. At higher activity one obtains a state of fully developed turbulence with a much higher density of defects (Figs. 2d and 3c). In this highly turbulent regime we expect the emergent length scale to become much smaller than, and therefore independent of, the system size . The above scaling form is then accordingly expected to reduce to
(16) 
In our simulations all scaling law measurements are taken safely within this regime of fully developed turbulence, such that the emergent length scales are free of finite size effects. We also explicitly demonstrate that finitesize effects indeed return when is no longer small, as illustrated by the snapshot of Fig. 2c.
It is also worth noting that at extremely large activities the defect density could in principle become so large that the defect spacing approaches the microscopic length scale . In this regime we would expect to be unable to decrease further upon any additional increase in activity, and so to saturate. However our simulations do not reach this limit and the inequality is always respected.
4.4 Form of the scaling function
Having proposed the existence of a scaling function in Eqn. 16, we now consider possible specific forms for this functional dependence of on the model parameters. Conflicting scaling laws for have been proposed in the existing literature ^{15, 16, 17}. While all of these studies agree that , there remains an apparent discrepancy over the scaling of with the activity.
Using Model II, Thampi et al. ^{15, 16} have proposed that , which would correspond to in Eqn. 16 having a square root dependence on its first argument. In contrast, using Model I, Giomi ^{17} suggested the relation , which would correspond to a linear dependence of on its first argument.
A possible origin of this discrepancy is the differing dimensionality of the order parameter between the two studies: while both have , Refs. ^{15, 16} had , whereas Ref. ^{17} had . This motivates us to compare numerical results for both within a single study. However our results below will rule out differences in as a source of discrepancy. Another potential reason could be that the two studies in fact explored different parameter regimes given the high dimensionality of the parameter space in these models. Therefore in order to ascertain the generality of these scaling laws, we systematically explore wide ranges for the three relevant adimensional parameters for both models. Our results will show that both forms suggested by the earlier studies can indeed apply, each in a different regime: one in the regime of fully developed active turbulence, the other when the system size plays a nontrivial role.
5 Results
We now present the results of our simulations. We focus on the regime of fully developed turbulence, corresponding to activity large enough to avoid finite systemsize effects () and yet small enough to avoid saturation of the defect spacing at the microscopic length (). We systematically explore the functional dependence of the emergent correlation lengths defined in Section 4.2 on the model parameters. Specifically in Model I we vary the activity, , and viscosity ratio, , keeping all other parameters fixed to the values in Table 2. In Model II we vary the activity and the nematic persistence length , with all other parameters fixed to the values in Table 2. We will show that in the region of fully developed active turbulence all of the emergent length scales defined above scale with the active length , in both models. We will additionally demonstrate that a weaker exponent might be obtained in the regime of less well developed turbulence, where the typical size of the emergent structures is an appreciable fraction of the box size.
5.1 Correlation lengths
In this section, we present our results for the correlation lengths defined in Sec. 4.2. Our main focus will be on an extensile nematic, corresponding to . We shall briefly discuss the contractile case at the end of this section.
Extensile active matter
Orientational correlations. We begin by considering correlations in the nematic order parameter . Figs. 2a and 3a shows the director correlation length and the defect spacing as obtained from Model II and Model I, respectively. For sufficiently large activity , we find that in both models both lengths obey a clear scaling law (black dashed lines). Note that the defect spacing correlation length is consistently larger than by a factor . This is to be expected as correlations at the halfway point between two defects () should be similar to those at .
At smaller activities (i.e., for ) the data obtained with Model II show a saturation in the power law (leftmost data points in Fig. 2a). This can be attributed to that fact that the length scale of nematic structure now spans an appreciable fraction of the system size, as seen in the snapshots of Fig. 2c. It is possible that fitting a power law in this saturation regime could result in a less negative exponent than the found in the regime of fully developed turbulence. Indeed we find that the scaling suggested by Thampi et al. (purple dashed dotted line in Fig. 2a) matches our data reasonably well in this regime.
The data in Fig. 2a also suggests that both and scale linearly with . We verify this scaling explicitly in Fig. 2b by plotting and against activity. The data for various values of collapse neatly onto a single curve, demonstrating a clear linear relation between both correlation lengths and .
The data obtained with Model I shown in Fig. 3a focus on large activities and verify that in this regime the scaling of both and with holds regardless of the model used. (They do not probe the saturation with system size seen at lower activities in Model II.) Data obtained for different values of the viscosity ratio can be collapsed when plotted as shown in Fig. 3b, suggesting , although the range of variation of the viscosity ratio is not sufficient to provide convincing evidence of scaling.
Taken together, the data obtained with the two models tests the functional dependence of the two nematic correlation lengths with respect to activity and the nematic persistence length . Once free of the system size, we find that both obey . Consistent with this scaling, replotting in Fig. 1 (main) the full director correlation function as a function of the rescaled coordinate gives good data collapse. Additionally, the data obtained with Model I suggest a scaling , but a larger range of values would be needed to verify this. Next we demonstrate that the same scaling form is observed for correlations lengths associated with the velocity field .
Velocity and vorticity correlation lengths. Using data obtained with Model I, we explore the dependence of the velocity correlation length and the vorticity correlation length (as defined in Section 4.2) on activity. In light of the results of the previous section, we directly plot both these lengths against the rescaled activity (see Fig. 4a). As shown previously for the orientation correlation lengths, we again observe that both and scale as , with all data sets falling approximately on a single curve. We stress that this behavior is different from that reported in Ref. ^{15}, where it was argued that does not depend on activity, while scale as .
Contractile active matter.
So far we have presented data for extensile systems, corresponding to . However many examples of contractile active matter are found in nature, e.g., suspensions of Chlamydomonas algae^{33}, or cytoskeletal actomyosin networks ^{34}. Therefore in order to further demonstrate the generality of our results, we now briefly consider the contractile case (). Since the linear instability of the homogeneous state requires requires ^{13}, for contractile systems we use . Our data, shown in Figs. 4(b,c), support the idea that the defect spacing (), director correlation length (), velocity () and vorticity () correlation lengths are all controlled by a single active length scale . We caution, however, that the mapping of rodlike extensile () onto disclike contractile () only holds at the linear instability level; the full nonlinear dynamics may be subject to additional instabilities depending on the specific parameter values. While we do not expect that this would significantly change the scaling behaviour, we defer a full study of these effects to future work.
5.2 Kinetic energy and enstrophy
The above scaling relations were obtained using the correlation functions defined in Section 4, which are normalised so that each function, e.g., , approaches unity as the separation distance . (See Fig. 1.) The normalization constants themselves, however, (i.e., the denominators in Eqs. 13 and 14) also provide useful information as they are directly proportional to the mean kinetic energy and enstrophy of the system, given by
(17)  
(18) 
where the angular brackets again denote an average over space and time. These quantities can be obtained experimentally, for instance by using particle image velocity (PIV) to quantify the flow fields of active liquids, as done by Dunkel et al. ^{35} in suspensions of extensile B. subtilis bacteria. We now use our earlier findings to motivate the expected scaling relation of these flow properties with activity, and then verify our predictions with further numerical data from both models.
Using simple dimensional analysis, the characteristic velocity of activityinduced shear flows associated with distortion of the local nematic order over a length scale can be obtained from the force balance condition (Eq. 1) as
(19) 
Our results indicate that for sufficiently large values of activity the physics is controlled by a single active length scale , with . Using this in Eq. 19, we find and . The scaling of the vorticity can be estimated as , which gives an enstrophy .
This scaling is consistent with the findings of Ref. ^{17} in which the author examined the typical size of vortex structures in the regime of spatiotemporal chaotic dynamics using what we refer to here as Model I and found that both the vortex size and the defect spacing appear to scale with the active lengthscale . Further evidence for this scaling can be found in the experiments of Ref. ^{35}, which found that where is the characteristic vortex size: assuming that , this implies that as we have argued above. Our proposed scaling is not, however, in agreement with the findings of Ref. ^{14, 16}. In those studies it was found that and , a result that cannot seemingly be reconciled with the simple assumption that .
In order to appraise these conflicting scaling laws, we perform simulations with both Model I and Model II and measure the kinetic energy (, see Figs. 5a, b) and enstrophy (, Figs. 5c, d). The data from both models clearly obey our expected scaling laws and (black dashed lines). With Model I our choice of units means that increasing the rotational viscosity is equivalent to reducing the solvent viscosity . Our dimensional analysis in Eq. 19 suggests that increasing should increase the characteristic velocity. We indeed observe this trend in our data in Figs. 4b, although simulations over a larger range of would be required to determine the exact scaling. By the same analysis, we also expect that should be proportional to for fixed , since . Our data from Model II explores several values of , and plotting against activity indeed leads to a reasonable curve collapse (Fig. 5d inset). Consistent with the findings of Giomi^{17}, we observe no appreciable dependence of on . This follows again from the scaling, .
6 Discussion
Using two distinct continuum models that have been studied extensively within the literature, we have performed a detailed numerical study of an active nematic to examine the scaling with activity of a number of structural and hydrodynamic correlation lengths, including the mean defect spacing. Our findings are consistent with the suggestion first put forward in Ref. ^{17} that in the regime of fully developed active turbulence defect proliferation, and the associated turbulentlike dynamics of the active nematic, are controlled by a single length scale . This is also the length scale that controls the onset of spontaneous flow instability of active films ^{2, 3, 36}. Our numerical data from both models show that all measures of correlation length considered scale with this length scale, for both extensile and contractile systems.
Two caveats must, however, be applied. First, for extremely large activities (i.e., ) activityinduced deformations below the nematic persistence length are expected to be suppressed. Secondly, at low activities, structures can form that span the system size, and correlation lengths will correspondingly saturate, (i.e., ). We have explicitly demonstrated this systemsize saturation in our simulations, a result that reconciles the apparently conflicting power law exponents previously reported in the literature.
Finally, to further support our findings, we have calculated the average kinetic energy and average enstrophy of the system, quantities that are readily obtainable from experiment. Our numerical results show that the scaling of these quantities with activity is consistent with a simple dimensional analysis based on the assumption that the physics is controlled by the single length scale .
Our results show that the key scaling relations hold for both strictly 2D and quasi2D models. Encouragingly, this implies that such models capture the dynamics of active nematics in a generic way, i.e., independent of the specifics of the model. How our results would compare with the equivalent fully simulation of an active nematic remains an interesting open question.
7 Acknowledgment
We thank Mark Bowick for introducing us to Ref. ^{31} (that describes a method for defect tracking) and Luca Giomi for writing the code used to study Model I. We thank them both and Mike Cates for invaluable discussions. MCM and PM were supported by the National Science Foundation through award DMR1305184 and by the Syracuse Soft Matter Program. EJH thanks EPSRC for a Studentship. SMF’s and EJH’s research leading to these results has received funding from the European Research Council under the European Union’s Seventh Framework Programme (FP7/20072013) / ERC grant agreement number 279365. The authors thank the KITP at the University of California, Santa Barbara, where they were supported through National Science Foundation Grant No. NSF PHY1125925.
Appendix A Flowaligning Parameter
We include here a comparison between the tumbling parameter used here and the LeslieErickson (LE) tumbling parameter (where corresponds to flowaligning regime and corresponds to flowtumbling regime). This comparison is presented in Appendix B of Ref. ^{3} for the case , but to our knowledge has not been displayed before for the case .
In dimensions the nematic tensor of a uniaxial nematic can be written as
(20) 
where and are the sum and difference, respectively of the two principal values of the moment of inertia tensor of uniaxial nematogens. In our case we use that corresponds to needlelike molecules^{37}.
We write the dynamical equation for the alignment tensor in dimension using the notation of Olmsted^{29},
(21) 
where , , and are parameters that couple order and flow.
Substituting the expression given in Eq. 20 for the alignment tensor into the dynamical equation Eq. 21, and assuming to be constant, we obtain an equation for the director,
(22) 
Comparing Eq. 22 to the LeslieErickson equation ^{29},
(23) 
we identify the correspondence between the Olmsted coefficients and the LeslieErickson coefficients as
(24) 
(25) 
Using in Eq. 24, we obtain
(26) 
Finally, for the case , we find
(27) 
The case was previously reported in Ref. ^{3}.
Footnotes
 footnotetext: Department of Physics, Durham University, Science Laboratories, South Road, Durham, DH1 3LE, UK
 footnotetext: Physics Department and Syracuse Soft Matter Program, Syracuse University, Syracuse, NY 13244, USA
 footnotetext: † Electronic Supplementary Information (ESI) available: videos of the two regimes of activity driven turbulence in Model II. See DOI: 10.1039/C6SM00812G.
 footnotetext: ‡ These authors contributed equally to this work.
 Note that the typical activityinduced velocity scale is in the range (see Fig. 5b or Eq. 19), meaning that the effective Reynolds number is in the range . These values are still small enough to ensure that any turbulence is activity driven (rather than inertial) in nature.
 In Model I the isotropicnematic transition is continuous and the equilibrium nematic correlation length given by diverges at the transition. Deep in the nematic state where we can approximate .
References
 M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao and R. A. Simha, Rev. Mod. Phys., 2013, 85, 1143–1189.
 R. Voituriez, J. F. Joanny and J. Prost, Europhys. Lett., 2005, 70, 404–410.
 D. Marenduzzo, E. Orlandini, M. E. Cates and J. M. Yeomans, Phys. Rev. E, 2007, 76, 031921.
 L. Giomi, M. C. Marchetti and T. B. Liverpool, Phys. Rev. Lett., 2008, 101, 198101.
 S. Ramaswamy, R. A. Simha and J. Toner, EPL (Europhysics Letters), 2003, 62, 196.
 S. Mishra and S. Ramaswamy, Phys. Rev. Lett., 2006, 97, 090602.
 V. Narayan, S. Ramaswamy and N. Menon, Science, 2007, 317, 105.
 L. Giomi, L. Mahadevan, B. Chakraborty and M. F. Hagan, Phys. Rev. Lett., 2011, 106, 218101.
 L. Giomi, L. Mahadevan, B. Chakraborty and M. F. Hagan, Nonlinearity, 2012, 25, 2245.
 L. M. Pismen, Phys. Rev. E, 2013, 88, 050502.
 L. Giomi, M. J. Bowick, X. Ma and M. C. Marchetti, Phys. Rev. Lett., 2013, 110, 228101.
 Shi Xiaqing and Ma Yuqiang, Nat Commun, 2013, 4, 1218.
 L. Giomi, M. J. Bowick, P. Mishra, R. Sknepnek and M. Cristina Marchetti, Philosophical Transactions of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 2014, 372, 0365.
 S. P. Thampi, R. Golestanian and J. M. Yeomans, Phys. Rev. Lett., 2013, 111, 118101.
 S. P. Thampi, R. Golestanian and J. M. Yeomans, EPL (Europhysics Letters), 2014, 105, 18001.
 S. P. Thampi, R. Golestanian and J. M. Yeomans, Philosophical Transactions of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 2014, 372, 0366.
 L. Giomi, Phys. Rev. X, 2015, 5, 031003.
 S. P. Thampi, R. Golestanian and J. M. Yeomans, Molecular Physics, 2015, 113, 2656–2665.
 T. Sanchez, D. T. N. Chen, S. J. DeCamp, M. Heymann and Z. Dogic, Nature, 2012, 491, 431–434.
 F. C. Keber, E. Loiseau, T. Sanchez, S. J. DeCamp, L. Giomi, M. J. Bowick, M. C. Marchetti, Z. Dogic and A. R. Bausch, Science, 2014, 345, 1135.
 G. Henkin, S. J. DeCamp, D. T. N. Chen, T. Sanchez and Z. Dogic, Philosophical Transactions of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 2014, 372, 0142.
 S. Zhou, A. Sokolov, O. D. Lavrentovich and I. S. Aranson, Proc. Nat. Acad. Sci. U.S.A., 2013, 111, 1265–1270.
 G. Duclos, S. Garcia, Y. HG and P. Silberzan, Soft Matter, 2014, 10, 2346–2353.
 T. Gao, R. Blackwell, M. A. Glaser, M. D. Betterton and M. J. Shelley, Phys. Rev. Lett., 2015, 114, 048101.
 P. de Gennes and J. Prost, The Physics of Liquid Crystals, Oxford: Oxford University Press., 2nd edn., 1993.
 A. Sokolov, I. S. Aranson, J. O. Kessler and R. E. Goldstein, Phys. Rev. Lett., 2007, 98, 1–4.
 H. H. Wensink, J. Dunkel, S. Heidenreich, K. Drescher, R. E. Goldstein, H. Lowen and J. M. Yeomans, Proc. Natl. Acad. Sci., 2012, 109, 14308–14313.
 M. Ravnik and J. M. Yeomans, Phys. Rev. Lett., 2013, 110, 026001.
 P. Olmsted, Rheologica Acta, 2008, 47, 283–300.
 E. J. Hemingway, A. Maitra, S. Banerjee, M. C. Marchetti, S. Ramaswamy, S. M. Fielding and M. E. Cates, Physical Review Letters, 2015, 114, 098302.
 D. Huterer and T. Vachaspati, Phys.Rev. D, 2005, 72, 043004.
 S. M. Fielding, D. Marenduzzo and M. E. Cates, Physical Review E, 2011, 83, 041910.
 S. Rafaï, L. Jibuti and P. Peyla, Phys. Rev. Lett., 2010, 104, 1–4.
 P. M. Bendix, G. H. Koenderink, D. Cuvelier, Z. Dogic, B. N. Koeleman, W. M. Brieher, C. M. Field, L. Mahadevan and D. A. Weitz, Biophys. J., 2008, 94, 3126–3136.
 J. Dunkel, S. Heidenreich, K. Drescher, H. H. Wensink, M. Bär and R. E. Goldstein, Phys. Rev. Lett., 2013, 110, 228102.
 S. A. Edwards and J. M. Yeomans, EPL (Europhysics Letters), 2009, 85, 18008.
 H. Stark and T. C. Lubensky, Phys. Rev. E, 2003, 67, 061709.