Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Transitions in the collective behavior of microswimmers confined to microfluidic chips
(USC Thesis Other)
Transitions in the collective behavior of microswimmers confined to microfluidic chips
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
Transitions in the Collective Behavior of Microswimmers Conned to Micro uidic Chips by Alan Cheng Hou Tsang A Dissertation Presented to the FACULTY OF THE GRADUATE SCHOOL UNIVERSITY OF SOUTHERN CALIFORNIA In Partial Fulllment of the Requirements of the Degree DOCTOR OF PHILOSOPHY (Aerospace and Mechanical Engineering) May 2016 Copyright 2016 Alan Cheng Hou Tsang Dedication To my parents ii Acknowledgements First of all, I would like to thank my thesis advisor, Eva Kanso, for her en- couragement, patience and trust. Her faith and optimism that I would be able to x the research problems have given me energy and enthusiasm to overcome the many research diculties I encountered throughout my doctoral study. It is always a great pleasure for me to work with Eva, and the things I have learnt from her are far more than just academic knowledge. Eva will remain my role model of scientist and teacher for the rest of my career. I want to thank my thesis committee members Paul Newton, Mitul Luhar, Christoph Haselwandter, and Peter Kuhn for their advice and guidance. I would also like to thank Larry Redekopp and Satwindar Sadhal for for being on my qualifying exam committee. Paul and Larry also taught several classes I enjoyed most at USC. I thank all my labmates for making my doctoral life so much fun. First, I thank Hanliang Guo for being a great ocemate since my rst semester at USC. Hanliang has been my very good friend and my best lunch buddy. I truly appreciate his support throughout these years. I thank Yangyang Huang who will get her PhD at the age I started my PhD life. Her achievement at such a young age always reminds me to work harder. I thank Sarine Babikian for always having funny stories to share. I thank Lionel Vincent for having grateful discussions with me on various topic of uid dynamics and his experiments have always been a source of amazement and inspiration to me. I thank Brendan Colvert for having interesting discussions with me on academia and education. I thank Anthony Medrano for being my rst iii mentee, and for asking questions that I don't have the answers to. I also thank the old members of the lab, Andrew Tchieu, Fangxu Jing, Babak Oskouei, Yang Ding, and Johannes Rudolph. I really enjoyed sharing the oce with them. I thank all my other friends at USC, Mak Tsang, Zaki Hasnain, Jerey West, Jeremy Mason, Angie Lee, Mohammad AlHamli, Dejuan Kong, Komsan Rattanaki- jsuntorn, Xinjiang Xiang, Trystan Madison, Prabu Sellappan. Special thanks go to Mohammad AlHamli for encouraging me to pursue an academic career. Most of all, I thank my parents for their love and patience. iv Table of Contents Dedication ii Acknowledgements iii List of Figures vii Abstract xiv Chapter 1 Scientic context and scope 1 1.1 Thesis overview . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1 1.2 Literature overview: experiments . . . . . . . . . . . . . . . . . . . 6 1.3 Literature overview: theories and simulations . . . . . . . . . . . . . 11 Chapter 2 Background 16 2.1 Hele-Shaw set-up . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.2 Balance of uid drag and wall friction . . . . . . . . . . . . . . . . . 19 2.3 Driven particle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 2.4 Swimming particle . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 Chapter 3 Problem formulation 28 3.1 Symmetric Microswimmer Model . . . . . . . . . . . . . . . . . . . 30 3.2 Asymmetric Microswimmer Model . . . . . . . . . . . . . . . . . . . 35 3.3 Potential dipole in conned systems . . . . . . . . . . . . . . . . . . 39 Chapter 4 Interaction of a pair of conned microswimmers 45 4.1 Synchronization, divergence and collision of symmetric swimmers . 46 4.2 Pursuit and synchronization of asymmetric swimmers . . . . . . . . 52 4.3 Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60 v Chapter 5 Flagella-induced transition in the collective behavior of conned microswimmers 62 5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 5.2 Modeling agellar activity . . . . . . . . . . . . . . . . . . . . . . . 64 5.3 Hydrodynamic model for conned microswimmers in doubly-periodic domain . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67 5.4 Emergence behaviors: chaotic swirling, global orientational order, and stable aggregation . . . . . . . . . . . . . . . . . . . . . . . . . 70 5.5 Statistics of global modes . . . . . . . . . . . . . . . . . . . . . . . 74 5.6 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79 Chapter 6 Circularly conned microswimmers exhibit multiple global patterns 81 6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 6.2 Hydrodynamic model for conned microswimmers in circular domain 83 6.3 Emergence behaviors: chaotic swirling, global vortex circulation, boundary aggregation . . . . . . . . . . . . . . . . . . . . . . . . . . 88 6.4 Statistics of global modes . . . . . . . . . . . . . . . . . . . . . . . 91 6.5 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 Chapter 7 Density shock waves in conned microswimmers 97 7.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97 7.2 Hydrodynamic model for conned microswimmers in rectangular chan- nel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99 7.3 Speed of sound . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 7.4 Subsonic to supersonic transition in density shock waves . . . . . . 108 7.5 Continuum equation model for density shock formation . . . . . . . 116 7.6 Control density distribution and population speed . . . . . . . . . . 118 7.7 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119 Chapter 8 Conclusions 122 References 127 Appendices 140 Appendix A Derivation of ow velocity eld of potential dipole in a doubly periodic domain 141 Appendix B Evaluation of Weierstrass type elliptic functions 145 vi List of Figures 1.1 Examples of emergent collective behaviors of active matter at dier- ent length scales: (a) self-organization of sh schools into a big vor- tex [135]; (b) ocking patterns formed by thousands of birds [135]; (c) microbial biolms in a sludge granule [29]; (d) motile bacteria exhibiting turbulentlike collective swimming [34]; (e) dynamic clus- ters of bacteria moving cooperatively [148]; (f) swirling pattern of coherently moving actin laments [108]; (g) emergent vortex lattice from moving microtubules [121]; (h) aggregation of light-activated colloids into crystal structures [92]. The length scales of the systems are indicated by the scale bars, above which the actual or approx- imate length of the systems are shown. The gures are adapted from the references cited above. . . . . . . . . . . . . . . . . . . . 3 1.2 Emergence of stable vortex state in circular connement: (a) motile bacteria self-organize into counter-rotating vortices, constituting a single vortex of swimmers encircled by a counter-rotating layer of swimmers [77]; (b) self-propelled colloidal rolling particles self- organize into a rotating vortex, characterized by a low density cen- ter and a high density vortex edge, as shown in the time-averaged local density eld [11]. The gures are adapted from the refer- ences cited above. . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 vii 1.3 Propagation of density sound waves in conned micro uidic droplets that are driven by external background ow. (a) Sound waves are observed in the form of `phonons' in a 1D micro uidic droplet crys- tal [7]. The normal modes of the phonons are given by the traveling transverse and longitudinal waves. The dispersion relations of the transverse and longitudinal modes are give by the peak of the in- tensity plot of the power spectrums in the bottom panels, where ! is the wave frequency and k is the wavenumber. (b) Sound waves propagate in a 2D ensemble of micro uidic droplet [32]. Disper- sion relation of the sound waves are depicted in the bottom panels, where ! and q x correspond to the wave frequency and the trans- verse wavenumber, and is the area fraction of the droplets in the channel. The gures are adapted from the references cited above. 9 1.4 Propagating density bands and density shock waves of conned active and driven particles in micro uidic channel. (a) A popula- tion of colloidal rolling particles exhibit unidirectional motion and emerge as a density band [12]. (b) Density shock evolves from a lo- calized jam of 1D micro uidic droplets driven by external ow [17]. The shock is located at the back of the droplets' jam. (b) Den- sity shock formed from a plug of driven 2D micro udic droplets [8]. The shock is formed at the front of the droplets. The gures are adapted from the references cited above. . . . . . . . . . . . . . . 10 1.5 Phase diagram of emergent collective modes for a population of asymmetric microswimmers conned in a Hele-Shaw cell. The re- sults are predicted based on the linear long-wave instabilities of the kinetic model introduced in [13]. The gures are adpated from [13]. 15 2.1 (a) Schematic of a Hele-Shaw cell. (b) Side view of the Hele-Shaw cell: the ow passing trough the conning plates has a Poiseuille prole. (c) Top view of the Hele-Shaw cell: the streamlines in any layers of the Hele-Shaw cell are identical and obey the potential ow. 17 2.2 Unconned and conned microswimmers have distinct hydrody- namic signature in far eld in the active and driven system. (a) In the active system, the far-eld in unconned and conned ow induced by swimmer's motion is given by force dipole and potential dipole respectively, (b) whereas in the driven system, the far-eld in unconned and conned ow induced by uniform background ow is given by stokeslet and potential dipole respectively. . . . . 23 viii 3.1 Dipole model for conned microswimmers: (a) symmetric (blu and slender) swimmer; (b) asymmetric (large tail and large head) swimmer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 3.2 (a) Flow velocity eld and (b) streamline of a potential dipole in an unbounded domain. The ow close to the singularity is covered. 40 3.3 (a) Flow velocity eld and (b) streamline of a potential dipole and its closest images in a doubly-periodic domain. The ow close to the singularity is covered. . . . . . . . . . . . . . . . . . . . . . . . 41 3.4 (a) Flow velocity eld and (b) streamline of a potential dipole in a conned circular domain. The circular boundary is denoted in red color. The ow close to the singularity is covered. . . . . . . . . . 42 3.5 (a) Flow velocity eld and (b) streamline of a potential dipole in a rectangular channel. The red line indicates the sidewall boundaries of the channel. The ow close to the singularity is covered. . . . . 44 4.1 Blu swimmer pairs exhibit (a) collision, (b) synchronization and (c) divergence of their paths. Parameter values: = 0:5, 1 =0:5, z 1 (0) = 1 (0) = 2 (0) = 0. (a)z 2 (0) = 5i; (b)z 2 (0) = 1 + 2:5i; and (c)z 2 (0) = 5 + 2:5i. Total integration time: (a) 50; (b) 60; (c) 100. The positions of the dipoles are marked by `' at t = 0 and by `o' at the end of the integration time. . . . . . . . . . . . . . . . . . . 47 4.2 Slender swimmer pairs exhibit (a) synchronization and (b) diver- gence of their paths. Parameter values: = 0:5, 1 = 0:5,z 1 (0) = 1 (0) = 2 (0) = 0. (a)z 2 (0) = 2+2i; (b)z 2 (0) = 2i. Total integra- tion time: (a) 90; (b) 30. The positions of the dipoles are marked by `' at t = 0 and by `o' at the end of the integration time. . . . 48 4.3 Poincar e section for (a) the case shown in Fig. 4.1(b) and (b) the case shown in Fig. 4.2(a). Note that the Poincar e sections (x 1 x 2 ;y 1 y 2 ; 1 ) and (x 1 x 2 ;y 1 y 2 ; 2 ) are coincident because 1 and 2 for these trajectories dier by a phase shift only. . . . . . . 48 4.4 Blu (top row) and slender (bottom row) swimmer: phase diagrams showing the behavior of the blu and slender swimmer pair as func- tion of initial conditions. The horizontal and vertical axis of the diagrams correspond to the real and imaginary part of the initial relative spacing of the swimmers z 2 (0)z 1 (0) respectively. The black regions correspond to initial conditions where the two swim- mers diverge, the grey regions correspond to synchronized motion of the two swimmers, and the white regions (except the excluded white disc around the origin) correspond to colliding trajectories. Note that x and y in the plots are in units of p . . . . . . . . . . 49 ix 4.5 (a) Two modes of solutions are obtained from (4.7): a pursuit mode (blue) where the swimmers trail one another and a synchronization mode (red) where they swim side by side. The shown solid curves correspond to ! 1 =i! 2 = 5 and the dashed lines to ! 1 =i! 2 = 10. The separation distance c is set to c = 4. (b) Summary of the stability analysis for these two modes. . . . . . . . . . . . . . . . . 53 4.6 Aperiodic behavior of two dipoles in doubly-periodic domain. The parameter values are ==3,z 1 (0) =2 exp(i),z 2 (0) = 2 exp(i), while is obtained by numerically solving the rst condition in (4.7). The positions of the dipoles are marked by `' at t = 0 and by `o' at the end of the integration. As time progresses, the two trajectories densely ll the whole domain. . . . . . . . . . . . . . . 55 4.7 (a) Pursuit and (b) synchronization in two dipoles undergoing pe- riodic motion. The parameter values are = tan 1 (3), z 1 (0) = 2 exp(i), z 2 (0) = 2 exp(i), while is obtained by numerically solving the rst condition in (4.7). The positions of the dipoles are marked by `' att = 0 and by `o' at the end of the integration time t = 80. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57 4.8 Pursuit and synchronization modes as attracting modes. Two dipoles hone in on quasi periodic trajectories where: (a) one dipole is in pursuit of the other for 2 = 0:5. (b) the two dipoles synchronize and move along parallel trajectories for 2 =0:5. The initial con- ditions are z 1 (0) = 1 (0) = 2 (0) = 0 while z 2 (0) = 1:5 + 1i in (a) and z 2 (0) = 2 + 1i in (b). . . . . . . . . . . . . . . . . . . . . . . . 58 5.1 Kinematic for the planar waving motion of a beating agellum in the laboratory frame and the swimming frame. Fundamental sin- gularities of potential ow, illustrated by the black dots, are dis- tributed along the agellum to calculate the ow induced by its beating motion. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65 5.2 (a) Time-averaged ow eld created by a beating agellum with A = 0:5, k = and U = 1. (Inset) Snapshots of the agellum- induced ow eld at dierent times. (b) Dipolar eld created by a circular disk of eective radius R tail tted to the average ow eld in (a). (c) Change in R tail with A and k of the traveling wave via the agellum. (d) Reduction of a agellated swimmer to a dumbbell swimmer. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 x 5.3 Strongly conned microswimmers in a Hele-Shaw cell create po- tential dipolar far-eld ows. Microswimmers with uniformly ran- domly oriented and spatially homogenous distribution in a doubly periodic domain are considered. . . . . . . . . . . . . . . . . . . . 68 5.4 (a) Swirling motion in vigorously-beating agellated swimmers; (b) orientational order; (c) aggregation in weakly-beating agellated swimmers. Probability distributions of velocity of the swimmers are depicted in the lower row. . . . . . . . . . . . . . . . . . . . . 71 5.5 (a) Swimmers with vigorously-beating agella orient with local ow and thus tend to chase each other; (b) Swimmers with weakly- beating agella orient in the opposite direction to local ow and thus tend to aggregate together. . . . . . . . . . . . . . . . . . . . 72 5.6 Statistics associated with swirling (top row) and orientational order (bottom row): polar order parameterhPi, velocity and angular cor- relations C v , C , probability distribution of rotational activity pa- rameter, all computed based on data corresponding to the shaded time range. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 75 5.7 (a) The time-averaged mean and standard deviation of and (b) the sample-averaged mean and standard deviation of long-time de- velopedhPi over various initial conditions. . . . . . . . . . . . . 77 5.8 Phase diagrams showing the three dierent collective modes, swirling, orientational order and aggregation, as a function of 1 , 2 and A . Regions of transitional behavior for which no clear pattern are iden- tied are left in white. . . . . . . . . . . . . . . . . . . . . . . . . 78 6.1 Microswimmers in Hele-Shaw and circular connement: (a) schematic of the set-up and the dipolar far-eld ows created by a swimmer in Hele-Shaw connement (inset); (b) Eect of the circular conne- ment is accounted for using an image system outside the circular domain as dictated by the Milne-Thompson circle theorem. . . . . 84 6.2 Snapshots for the emergence of three global modes: (a) chaotic swirling ( = 3); (b) stable circulation ( = 0:5); (c) boundary aggregation ( =0:5). The corresponding ow elds are depicted at the bottom row. For all three simulations, N = 750 and R = 50. 88 6.3 Statistical data for the stable circulation mode shown in Fig. 6.2(b): (a) vortex order parameter ; (b) mean tangential speed and prob- ability distribution of the swimmers; (c) mean polarization eld. The mean polarization eld and tangential speed are computed based on the averaging of data over the shaded time range, which a stable vortex order has been developed. . . . . . . . . . . . . . 93 xi 6.4 (a) The time-averaged mean and standard deviation of long-time developed over various initial conditions. The standard devi- ations are denoted by the error bars. (b) Phase diagram showing the three dierent global modes: chaotic swirling, stable circulation and boundary aggregation, as a function of and R. . . . . . . . 94 7.1 (a) Schematic of microswimmers in a Hele-Shaw rectangular chan- nel. (b) Dipolar far-eld ows induced by the self-propelled motion and the background ow. (c) Sidewall connement is accounted for using an innite image system. . . . . . . . . . . . . . . . . . . . 99 7.2 The gradient in the density distribution of the initial swimmer pop- ulation is essential for the formation of the compression shock waves.104 7.3 Speed of sound in the moving frame of reference with velocityUp+ V . Parameters used: d = 1, = 0:5. . . . . . . . . . . . . . . . 108 7.4 Density shock waves in conned microswimmers. Left: (a) for V < c, compression shock wave forms at the back; (b) for V = c, shock is suppressed; (c) for V >c, compression wave forms at the front. The `speed of sound' of the medium ahead is c = U +V . Right: Comparison of simulation results between discrete parti- cle and continuum model. The density is the local area fraction. Simulation results of the discrete particle model are colored in blue. Simulation results of the quasi 1D continuum model are superim- posed in red. Parameter values are: W = 60, A = 0:3, = 1, snapshots taken at t = 1000. . . . . . . . . . . . . . . . . . . . . . 110 7.5 Comparison of the speed of the density shock wave and the speed of sound for the case shown in Fig. 7.4. The probability distribution of swimming speeds of the microswimmers in the discrete simulations are shown in solid lines, and the speed of sound of the medium ahead is indicated in dashed lines. . . . . . . . . . . . . . . . . . 110 7.6 1D model does not capture correctly the shock wave formation. It results in compression shock at the front forV <c and at the back forV >c, as opposed to the 2D system which exhibits compression at the back for V <c and front for V >c (Fig. 7.4). . . . . . . . . 111 xii 7.7 (a) Shock order parameter is dened by comparing the actual den- sity prole to a triangular density prole, shown here att = 2000 for W = 80 and A = 0:3. (b) Evolution of the shock order parameter for the case in (a). (c) Mean value of the shock order parame- ter as a function of channel width W and swimmers' area fraction A , computed base on the shaded time range and averaged over all sample simulations. The dashed line shows the contour line of hSi t!1 = 0:4, above which shock waves are unambiguously observed.113 7.8 Mean value of the shock order parameter showing the transition of density shock waves from subsonic to supersonic via increasing background ow velocity V, averaged over several sample simula- tions. Parameters used: W = 60, = 0:5, A = 0:3, = 1. . . . . 115 7.9 The shock front of microswimmers propagates in a periodic channel lled with swimmers eventually, and evolve as a traveling density band. Parameter values are: W = 40, L = 400, A =0.3, = 0:5, = 1, V = 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 116 7.10 Microswimmers homogeneously lling the ow channel: a large den- sity perturbation propagates in the homogenous medium as a trav- eling density band. Parameter values are: W = 40, L = 400, A =0.3, = 0:5, = 1, V = 1:8, whereas A = 0:2 for the homogenous background density and A = 0:5 for the localized density perturbation. . . . . . . . . . . . . . . . . . . . . . . . . . 116 7.11 Segregation of two species of microswimmers that are mixed homo- geneously, with initial conditions given in (a). The two species of swimmers are of distinct (b) translational and (c) rotational motil- ity, respectively. Parameter values are: W = 80,L = 800, A = 0:2 for both species of swimmers. . . . . . . . . . . . . . . . . . . . . 120 7.12 The strength of the background ow is controlled to obtain a Gaus- sian density distribution of microswimmers from an initially uni- form density distribution. Parameter Values are: W = 80,L = 800, A = 0:3, = 0:5, = 1. Snapshots of the swimmers at dierent time are shown on the same axis for illustration. . . . . . . . . . . 120 xiii Abstract Micro uidics are revolutionizing modern research in disease diagnostics and drug development. For example, micro uidic devices for cell manipulation and sorting provide high-throughput methods for detection of infectious bacterial cells and oer great potential in cancer cell diagnostics. Great progress has been achieved in these directions, however, many questions concerning the design and control of these sys- tems remain open. In particular, there is limited understanding of the collective behavior of motile cells in conned micro uidic environments, which makes the pro- cess of controlling these cells particularly challenging. In this thesis, we investigate the collective behavior of cell populations in the conned micro uidic geometries. We develop physics-based, mathematical and computational models of bacterial populations in conned micro uidic environment. We particularly focus on the hydrodynamic interactions among members of these bacterial populations. As such, our models account for the long-range hydrodynamic signature of bacterial cells under connement and the short-range steric interactions among swimming bacteria and the micro uidic connement. In one set of models, we investigate the eects of agellar activity on the hy- drodynamic interactions in bacterial populations in a doubly-periodic domain. We show that decreasing agellar activity induces a hydrodynamically triggered tran- sition in conned bacteria from swirling to aggregation. The swirling behavior is consistent with recent experimental ndings, but our work goes beyond these nd- ings to provide a predictive tool for systematically analyzing the collective bacterial xiv behavior under various environmental conditions. In particular, our results show that varying the level of agellar activity, which can be accomplished by varying en- vironmental conditions, provide physical mechanisms for coordinating the collective behaviors in conned bacteria. This can in turn alter molecular microbial processes such as transport of oxygen and nutrients that mediate transitions in the bacteria physiological state. We then examine the eects of circular connement on these transitions of col- lective behavior. We show as in the case of doubly-periodic domain that decreasing agellar activity induces hydrodynamically triggered transitions, but from swirling to global circulation (vortex) to boundary aggregation and clustering. The spon- taneous self-organization of bacterial populations into global vortex is reminiscent to recent experimental observations on circularly-conned populations of motile bacteria and rolling colloids. These results highlight that the geometry of the envi- ronments have signicant in uences on the emergent global patterns. Finally, we show that bacterial populations subject to increasing external ow in micro uidic channels exhibit an interesting transition of density shock waves from `subsonic' to `supersonic' modes, characterized by the shift in location of high popu- lation density from the back to the front of the group. These ndings are consistent with experimental observations on driven droplets and self-propelled colloids. They also serve to guide the development of mechanisms for controlling the emergent den- sity distribution and the average population speed of motile cells, with potentially profound implications on the transport and sorting of cells in micro uidic channels. xv Chapter 1 Scientic context and scope 1.1 Thesis overview Active systems are systems driven internally by large collections of self-propelled individual units. They are ubiquitous in nature, with examples spanning a wide range of length scales, from schools of sh [25, 26] to populations of motile bacte- ria [21, 34, 36, 81, 114] and assemblages of sub-cellular extracts such as cytoskeletal networks [63, 64, 143, 144], actin laments [108] and microtubules [58, 87, 107, 121, 123]. These active systems often exhibit rich collective behavior at a system scale that is typically several orders of magnitude larger than the scale of the individual unit. In active systems, each individual unit usually moves according to its own self-propelled motion and in response to signals from the environment. The self- propelled motion is the result of energy input into the system by the individual unit. It usually comes from the storage of each individual unit or from an en- vironmental supply. The continuous injection of energy from the individual units 1 drives the system out of equilibrium, leading to the many surprising and fascinating large-scale collective behaviors through complex interactions among the individual units. A gallery of collective behaviors at dierent length scales in nature is shown in Fig. 1.1. Active systems are attractive model systems for scientic research and practical applications. In biological systems, collective behavior is believed to play impor- tant roles in the functioning and survival of the groups [135]. Inanimate systems such as driven and self-propelled droplets and reactive colloids provide an attractive paradigm for recongurable smart materials [85, 98, 122] and biomedical applica- tions such as cell testing and sorting in high-throughput micro uidic `lab-on-chip' devices [54, 84]. Further, the richness of dynamical patterns in active systems oer excellent examples to investigate fundamental physics problems such as non- equilibrium statistical mechanics, pattern formation, molecular diusion and active microrheology [15, 79, 97, 105, 128]. In this thesis, we focus on the population of microswimmers in viscously domi- nant ow at low Reynolds number. The question of how the highly-coordinated collective motions arises from piecewise interactions among the microswimmers has been the subject of intense research in the past few years. However, most of the theoretical and experimental studies have focused on the motion in three- dimensional (3D) uid environments where the eect of geometric connement to the microswimmers at a population level is weak or absent. Motivated by recent technological advances in producing and manipulating large ensembles of particles 2 (e) (d) (f) (g) (h) 10μm 25μm 10μm 50μm 2mm (a) (b) (c) ~ O (1m) ~ O (10m) 20μm Figure 1.1: Examples of emergent collective behaviors of active matter at dierent length scales: (a) self-organization of sh schools into a big vortex [135]; (b) ock- ing patterns formed by thousands of birds [135]; (c) microbial biolms in a sludge granule [29]; (d) motile bacteria exhibiting turbulentlike collective swimming [34]; (e) dynamic clusters of bacteria moving cooperatively [148]; (f) swirling pattern of coherently moving actin laments [108]; (g) emergent vortex lattice from moving microtubules [121]; (h) aggregation of light-activated colloids into crystal struc- tures [92]. The length scales of the systems are indicated by the scale bars, above which the actual or approximate length of the systems are shown. The gures are adapted from the references cited above. 3 in micro uidic devices, attention began to shift to the collection behavior of the mi- croswimmers conned in quasi-2D geometries. Understanding the role of geometric connement on the collective behavior of microswimmers is directly relevant for the design and control of micro uidic devices for cell manipulation and sorting. We develop physics-based mathematical and computational models of microswimmers conned in Hele-Shaw type geometries, that is, microswimmers strongly conned by upper and lower plates, such that their motion is restricted to a two-dimensional (2D) plane. We consider the hydrodynamic interactions in the populations of such conned microswimmers and examine the transitions in their collective behavior. We show that the conned microswimmers exhibit rich emergent behavior, includ- ing chaotic swirling, global order, stable aggregation, and density shock waves. The outline of this thesis is given as follow. For the rest of Chapter 1, we brie y review the literature concerning the ex- periments, theories and simulations of population of microswimmers, and highlight several important results in conned systems. Chapter 2 provides several important background concepts. The idea of Hele- Shaw ow is brie y reviewed. We then revisit the physical principles for the changes in mobility properties and far-eld hydrodynamic signature of the microswimmer due to the Hele-Shaw connement. Chapter 3 provides the problem formulation for interacting conned microswim- mers. We discuss in detail how the geometry of the microswimmer aects the orien- tation dynamics of the microswimmer and group the microswimmer models into two 4 categories, that is, symmetric and asymmetric microswimmer models. The far-eld ows induced by the motion of a microswimmer in dierent conned geometries are also presented. Chapter 4 focuses on the dynamics of a pair of symmetric and asymmetric microswimmers and identies their interaction modes. Chapter 5 considers the collective behavior of conned microswimmers in a doubly-periodic domain. We generalize the asymmetric microswimmer model pre- sented in Chapter 3 by relating the model with the eects of agellar activity of the microswimmers. We explore the emergent collective modes of the microswimmers and their statistical properties. Chapter 6 considers a micro uidic system composed of both Hele-Shaw and cir- cular connement, and shows that populations of microswimmers subject to such a connement can exhibit a rich variety of global patterns, including a global cir- culating vortex state which is reminiscent to those observed in recent experimental studies. Chapter 7 focuses on the formation of density shock waves in the populations of conned microswimmers driven by external background ow in a rectangular channel. We contrast the dierence in the hydrodynamic interactions in 1D and 2D system and reveal the physical mechanisms for the shock formation. Finally, in Chapter 8, we comment on the signicance of these ndings and provide a brief outlook for future research. 5 1.2 Literature overview: experiments The physics of natural and synthetic microswimmers is an active area of research both at the individual and collective levels over the last decades. At the individ- ual level, a great deal of attention is given to deciphering or devising mechanisms for self-propulsion of microswimmers [68, 91]. When a pair of microswimmers come together, they sometimes synchronize with each other and move coopera- tively [95, 124, 145]. At the population level, these microswimmers exhibit rich emergent dynamics at the system's scale, which makes striking dierences between the dynamics of an isolated swimmer and the dynamics of a swimmer as part of the population. Populations of motile bacteria such as Escherichia coli or Bacillus subtilis are observed to exhibit transient, recurring states of collective motion in the form of large-scale jets and swirling vortices in many experiments [21, 22, 28, 34, 36, 67, 81, 114, 115, 139], see for example, Fig. 1.1(d). This chaotic collective motion of bacte- rial population is usually referred as `bacterial turbulence', albeit the ow is actually at low Reynolds number. In the state of bacterial turbulence, the bacteria exhibit local, but not global orientational order with the others, and their motions result in large-scale vortical uid ows that lead to increase in speed of the individuals and enhanced diusion of passive particles in the bacterial bath [59, 65, 71, 146]. In ad- dition to bacterial turbulence, motile bacteria also exhibit other collective behaviors such as swarming [19, 93], dynamic clustering [148, 149], and biolms [29, 37, 82]. Moreover, the ows induced by a population of swimming bacteria were shown to 6 have signicant impact on the eective viscosity of the uid [43, 96, 113], and most recently, it was observed that the eective shear viscosity of the uid vanishes in a bacterial bath with semi-dilute concentration, thereby turning the viscous uid into a super uid [74, 78]. In term of practical applications, the self-organization of populations of bacteria can be utilized to develop bacterial-driven micro uidic devices, including micro uidic mixing systems, micro uidic pump and ratchet mo- tors [60, 61, 70]. Advances in nanotechnology and micro uidics has lead to great progress on analyzing the emergent behavior of large populations of motile colloids. These motile colloids are not autonomous, which means that they require external en- vironmental conditions such as chemical or surface tension gradient, magnetic or electric elds to deliver their motions. Therefore, the rule of interactions among the motile colloids vary in dierent systems. Populations of motile colloids ex- hibit a rich variety of collective behaviors. For example, motile rolling colloids and self-propelled droplets are observed to form ordered and directional motions spontaneously [11, 12, 126]. Reactive colloids are observed to self-assemble into crystal structures [14, 92, 110, 125, 138], providing new methods for synthesis of nanomaterials. Recent experiments demonstrated that geometric connement has striking ef- fects on the self-organization and emergence of global patterns for populations of motile colloids and driven micro uidic droplets. Here, we highlights the results of these experiments which motivate the studies in this thesis, including spontaneous 7 (a) (b) Figure 1.2: Emergence of stable vortex state in circular connement: (a) motile bacteria self-organize into counter-rotating vortices, constituting a single vortex of swimmers encircled by a counter-rotating layer of swimmers [77]; (b) self-propelled colloidal rolling particles self-organize into a rotating vortex, characterized by a low density center and a high density vortex edge, as shown in the time-averaged local density eld [11]. The gures are adapted from the references cited above. circulation of bacteria and motile colloids in circular connement [11, 77, 142], prop- agation of density sound waves in micro uidic droplet connement [5, 7, 32], and formation of density bands and density shock waves of motile colloids and driven droplets in a micro uidic channels [8, 12, 17]. Experiments on conned bacterial suspensions within a cylindrical droplet showed that if the radius of circular connement is below certain critical value, a homoge- nous distribution of bacteria with random orientation will spontaneously organize into a spiral vortex encircled by a counter-rotating cell boundary layer [77, 142], see Fig. 1.2(a). This global ordering of bacteria into vortex state is the result of the complex interplay between radial connement, hydrodynamic and steric inter- actions, evidenced by the simulations from a computational model of interacting microswimmers that reproduce the counter-rotating vortices observed in the exper- iments [77]. In another experiment on electrically powered self-propelled rolling 8 (a) (b) Figure 1.3: Propagation of density sound waves in conned micro uidic droplets that are driven by external background ow. (a) Sound waves are observed in the form of `phonons' in a 1D micro uidic droplet crystal [7]. The normal modes of the phonons are given by the traveling transverse and longitudinal waves. The dispersion relations of the transverse and longitudinal modes are give by the peak of the intensity plot of the power spectrums in the bottom panels, where! is the wave frequency and k is the wavenumber. (b) Sound waves propagate in a 2D ensemble of micro uidic droplet [32]. Dispersion relation of the sound waves are depicted in the bottom panels, where ! and q x correspond to the wave frequency and the transverse wavenumber, and is the area fraction of the droplets in the channel. The gures are adapted from the references cited above. colloids conned in a quasi-2D cylinder, it was shown that the colloids can also spontaneously self-organize into a single vortex [11], see Fig. 1.2(b). However, in contrast to the counter-rotating vortices in bacteria, the vortex state formed by the colloids is characterized by a single vortex with a low density center and a high density edge at the circular boundary. Motile and driven particles conned in a quasi-2D micro uidic channels exhibit interesting emergent behavior, from propagating density sound waves, density bands to density shock waves, see Fig. 1.3 and Fig. 1.4. The formation of these density 9 (a) (b) (c) Figure 1.4: Propagating density bands and density shock waves of conned active and driven particles in micro uidic channel. (a) A population of colloidal rolling particles exhibit unidirectional motion and emerge as a density band [12]. (b) Density shock evolves from a localized jam of 1D micro uidic droplets driven by external ow [17]. The shock is located at the back of the droplets' jam. (b) Density shock formed from a plug of driven 2D micro udic droplets [8]. The shock is formed at the front of the droplets. The gures are adapted from the references cited above. waves is attributed to the long-range dipolar interaction among the particles. In a one-dimensional (1D) micro uidic droplet crystal, the sound waves propagate in the form of traveling `phonons' which include both transverse and longitudinal nor- mal modes, see Fig. 1.3(a). Sound waves are also observed in a 2D homogenous suspension of droplets, see Fig. 1.3(b). These sound waves have well dened disper- sion relations that characterize their propagation properties, such as the phase and group velocities. A more remarkable emergent behavior is observed when external ow is applied to a dense plug of micro uidic droplets, that is, the formation of 10 density shock waves. Interestingly, the locations of the compression shocks are dif- ferent in the 1D and 2D system. In 1D system, the compression shock is formed at the back of droplets' jam, see Fig. 1.4(a), whereas the compression shock is formed at the front in the case of the 2D micro uidic droplets, see Fig. 1.4(b). This con- tradiction reveals the complexity of the mechanism responsible for the formation of density shock waves and the dierence in nature between 1D and 2D hydrodynamic interaction. 1.3 Literature overview: theories and simulations Many mathematical and computational models have been proposed to simulate the emergent collective behavior observed in nature. These models can be broadly distinguished into rule-based and physics-based. The rst rule-based model was proposed by Vicsek et al. in [134] to investigate the emergence of ordered motion in systems of particles that mimic biological rules of interaction. Many similar models were proposed afterwards, see [16, 18, 25, 26, 135] and the references therein. These models consider homogeneous individual units interacting locally based on rules of repulsion, alignment, and attraction to other individuals, and are capable of exhibiting realistic dynamics similar to those observed in nature, yet they do not elucidate the underlying physical mechanisms governing the emergent behaviors. Physics-based models for microswimmers are ubiquitous. One major challenge in developing such physics-based models for populations of microswimmers is to cap- ture the hydrodynamic interaction among the microswimmers. Classical modeling 11 and simulation methods of solid- uid interactions are computationally prohibitive for this purpose, as the problems of collective behavior involve thousands of inter- acting swimming bodies and long-time dynamical simulations. Alternatively, many reduced models based on the long-range hydrodynamic and short-range steric inter- actions were developed. Namely, when the swimmers are in motion, they create ow disturbances that aect the motions of other swimmers in far-eld hydrodynami- cally. When the swimmers are in close proximity, they employ collision avoidance mechanisms in the form of steric repulsion and sometimes also reorientation. These models are usually developed using either discrete particle or continuum approach. The discrete particle models take into account the salient features of the swim- mer's geometry, without including all its details. Examples include minimal dumb- bell swimmers constituted of two beads linked by a frictionless rod [24, 44, 45], swimmers with the shape of a slender rod [76, 77, 101, 104], spherical squriming swimmers [31, 38, 49, 50, 51, 52, 53, 72], and rotor particles [100, 117, 147]. In these models, each individual swimmer is assumed to self-propel without taking into account the details of the propulsive mechanism of the swimmer. The eects of propulsive mechanisms employed by the swimmers are re ected on the ows in- duced by the individual swimmers. According to their propulsive mechanisms, the swimmers are usually categorized into pusher or puller [102, 103]. In the case of bacterium, a pusher corresponds to swimmer that exerts a propulsive force near the rear of the its body, such as Escherichia coli and Bacillus subtilis, whereas a puller, such as Chlamydomonas reinhardtii, exerts a thrust at its front. The ow 12 induced by pusher or puller in far-eld is characterized by a force dipole singularity. The ow direction of the force dipole induced by the pusher and puller are exactly opposite to each other. These discrete swimmer models successfully capture the emergence of large scale ows and collective coherent structures observed in many experiments on bacterial suspensions [21, 34, 36, 81, 114]. In the continuum approach, instead of considering the population of microswim- mers as individual interacting particles, the conguration of the swimmer is modeled as a time varying probability density function, with the location and orientation of the swimmers as the dependent variables [105, 106]. The density function is coupled with the governing equation of the uid ow via a density dependent stress term exerted by the swimmers on the local uid, and the density function evolves with time according to the self-propelled motion and the local ow in return. The contin- uum equations are developed based on phenomenological models [3, 35, 36, 83, 112] or mean-eld kinetic models that accounts for the conservation of particle proba- bility distribution using Smoluchowski equation [47, 102, 103, 120]. Extension of these models were proposed to study other physical eects in addition to hydrody- namic interactions, such as external ow [90], chemotaxis [39, 62, 75], steric inter- actions [2, 41], air-liquid interface [80], and connement [40, 143]. One of the merits of the continuum approach is that it allows analytical insight into the macroscopic behavior of the swimmers via simplication of the equation models. Through linear stability analysis on the continuum models, it was shown that a homogenous pop- ulation of pushers is always unstable, whereas a homogenous population of pullers 13 is always stable [102, 103]. As a result, large spatiotemporal uctuations are ob- served in the population of pushers, which lead to the emergence of turbulentlike motion. This agrees with the experimental results that only pusher-type bacterial populations are observed to exhibit bacterial turbulence. Recently, there are growing number of theoretical models proposed to under- stand the collective behavior of microswimmers conned in quasi-2D geometries such as Hele-Shaw cell. Geometric connement changes drastically the nature of the hydrodynamic interactions among particles [6, 13, 27, 33, 111]. The long-ranged hydrodynamic interactions driven by the 3D force dipole are screened by the solid conning walls, making it subdominant in comparison with the 2D potential dipole arising from incompressibility. As a result, the long-range interactions among the conned swimmers can be obtained from the superposition of potential dipole sin- gularities. Moreover, the far-eld hydrodynamic signature induced by a uniform external ow on a conned microswimmer is also given by potential dipole [33], but the direction of this induced dipole is always opposite to the background ow, irrespective to the orientation of the swimmer. More remarkably, in [13], Brotto et al. proposed a head-tail asymmetric dumbbell microswimmer model, and showed that due to friction with the nearby walls, such swimmers reorient in response to both the local ow eld and its gradient, as opposed to the unbounded swimmers that only reorient in response to the local ow gradient. Further, they performed a linear stability analysis on a kinetic continuum model and predicted that the con- ned asymmetric swimmers can exhibit new collective modes that are not observed 14 Figure 1.5: Phase diagram of emergent collective modes for a population of asym- metric microswimmers conned in a Hele-Shaw cell. The results are predicted based on the linear long-wave instabilities of the kinetic model introduced in [13]. The gures are adpated from [13]. in unbounded systems, see Fig. 1.5. Motivated by this result, discrete particle sim- ulations on a desingularized potential dipole model were performed in [69] and the simulations showed the emergence of new collective modes, including polarization waves, active lanes, and trac jams. They also proposed a continuum trac ow model to explain the formation of trac jams in their system. In a recent study on populations of spherical squirming swimmers in quasi-2D connement, inter- esting phase separation and clustering behavior were observed [151]. In summary, microswimmers conned in Hele-Shaw geometries have distinct hydrodynamic sig- nature, orientation dynamics, and collective behavior compared to their unconned counterpart. 15 Chapter 2 Background This thesis is concerned with the collective dynamics of microswimmers that are con- ned in Hele-Shaw type geometries. Here, we revisit the idea of Hele-Shaw ow and consider the role of uid drag and wall friction on force balance of the microswim- mer in a Hele-Shaw cell. In particular, we derive the expression for the translational mobility coecient of the microswimmer. We then contrast the dierences between far-eld ow singularities induced by the mircoswimmer in unconned and conned systems, and show that potential source dipole is the dominant ow singularity in the conned systems. 2.1 Hele-Shaw set-up We brie y review the concept of the uid ow in a Hele-Shaw cell. A Hele-Shaw cell is a micro uidic device which consists of two parallel plates separated by a small and uniform gap with depth of h in the z direction and much larger dimensions in the x and y directions, see Fig. 2.1(a). Note that if the conning surfaces are not 16 h (a) (b) (c) x z y Figure 2.1: (a) Schematic of a Hele-Shaw cell. (b) Side view of the Hele-Shaw cell: the ow passing trough the conning plates has a Poiseuille prole. (c) Top view of the Hele-Shaw cell: the streamlines in any layers of the Hele-Shaw cell are identical and obey the potential ow. smooth, the eect of surface roughness can be accounted for locally by considering h as a spatially varying function, that is, h = h(x;y). Consider a Newtonian and incompressible ow at zero Reynolds number. The ow velocity eld u(x;y;z) is governed by the Stokes equations r u = 0; r 2 u =rp; (2.1) where is the dynamic viscosity and p is the pressure. In a Hele-Shaw cell, the gap between the the parallel plates are small compared to the dimensions in thex-y plane,h! 0, and this geometric constraint has imposed a much larger ow gradient in the z direction, namely, @=@z @=@x;@=@y for the ow velocity. Given this 17 constraint, it is reasonable to assume that the pressure depends onx andy only, that is,p =p(x;y) and therefore we can separate the variables in u(x;y;z), expressing it as the product of az-dependent function and a (x;y)-dependent function. Using the no-slip boundary conditions at z = 0 and z =h, we can show that the momentum equation in (2.1) reduces to u = z(hz) 2 rp: (2.2) Recall that potential ow is dened by a ow proportional to the gradient of a scalar velocity potential function. From (2.2), we can see that at any given z, the ow velocity in thex-y plane of a Hele-Shaw cell corresponds to a two-dimensional (2D) potential ow in which the velocity potential is proportional to p. The ow between the conning plates is characterized by a parabolic Poiseuille prole and zero vorticity in the z-direction, see Fig. 2.1(b). If we consider dierent layers of uid in thez-direction, the streamlines of the uid ow in any layers have the same form, albeit having dierent magnitude of the ow velocity, see Fig. 2.1(c) for the streamlines passing through a xed circular cylinder in a Hele-Shaw cell. Therefore, the ow can be treated as a 2D ow by considering the gap-averaged ow velocity hui z = Z h 0 udz = h 2 12 rp: (2.3) 18 By virtue of (2.2) and (2.3), we obtain from the continuum equation in (2.1) the eld equation for the ow in a Hele-Shaw cell, given by the Laplace's equation r 2 = 0; (2.4) where the velocity potential (x;y) =h 2 =(12)p(x;y). To conclude, the ow in a Hele-Shaw cell is governed by the Laplace's equation given by (2.4), which can be solved using potential ow theory. The Laplace's equation is complemented by impenetrability conditions if any physical boundaries exist in the x-y plane. 2.2 Balance of uid drag and wall friction Consider a microswimmer located at r o translating with a velocity _ r o in a Hele-Shaw cell. The size of the swimmer is of the same order as the height of the Hele-Shaw cell h such that the swimmer is strongly conned by the upper and lower walls, and the swimmer is therefore restricted to move on a two-dimensional plane. A swimmer subject to a gap-averaged local ow velocity u(r o ) is not advected at the same velocity as the local ow, but at a velocity of u(r o ), where is dened as the translational mobility coecient. This means that the swimmer actually \senses" a smaller ow velocity than the actual ow velocity. This screening of ow is based on the fact that in a Hele-Shaw cell, the swimmer is subject not only to viscous drag from the uid, but also to lubricated friction from the upper and 19 lower conning walls. Here, we derive the expression of the translational mobility coecient in terms of the drag and friction coecient of the swimmer; for more details, see also [6, 133]. We start by considering the forces acting on the swimmer. For simplicity, we consider the case when the swimmer has the shape of a circular disc with radius of R d , but the general result obtained from the analysis holds for swimmer with arbi- trary shape. A self-propelled swimmer is subject to a propelling force F s generated by its swimming motion. When the swimmer experiences a non-zero local ow, the pressure gradient sustaining the local ow u(r o ) exerts a force F h on the swimmer, which can be obtained from (2.3) and is thus given by F h =u(r o ); (2.5) where = 12R 2 d =h and is the dynamic viscosity of the uid. When the swimmer is in motion, it experiences a viscous resistive drag F d from the uid, F d =(_ r o u(r o )); (2.6) where is the drag coecient which depends on the size and shape of the swimmer. For a circular swimmer, the drag coecient has the form [133] = 4 h " 3R 2 d h 2 + p 12R d K 1 ( p 12R d =h) hK 0 ( p 12R d =h) # ; (2.7) 20 whereK 0 andK 1 are the modied Bessel functions of the second kind. From (2.7), it can be seen that is an increasing function of R d . The friction F f acting on the swimmer from the conning plates is F f = _ r o ; (2.8) where is the friction coecient. If we assume simple shear ow between the swimmer and the lubrication gaps, then = 2R 2 d =h gap ; (2.9) where h gap is the depth of the gap between the swimmer and the conning plates. Now, we consider the force balance on the swimmer, in the absence of inertia, the sum of the forces acting on the swimmer is zero, F s + F h + F d + F f = 0: (2.10) Substitute (2.5) { (2.8) into (2.10), and rearrange, we obtain, _ r o = U +u(r o ); (2.11) where U = F s =( + ) is the self-propelling velocity of the swimmer. Note that the dierence in force balance between a self-propelled swimmer and a driven particle is only the propulsive force. That is, all formulation remains the same except that 21 F s = U = 0 for the case of driven particle. The translational mobility coecient is given by = + + ; (2.12) where > due to the small lubrication gap,h gap h, and therefore< 1. From (2.7) and (2.9), it is straightforward to show that= is a decreasing function of the swimmer size R d , hence decreases with the increase in swimmer size. This result agrees with the intuition that it needs more eort and hence larger ow to move a larger body between conning plates with the same amount of distance moved by a smaller body. We conclude this section by noting that the force free-equation given by (2.11) can be generalized to the case of N swimmers. This leads to the equations governing the translational motions for a population of microswimmers. For the n-th swimmer, n = 1;:::;N, the translational equation is given by _ r n = U +u(r n ): (2.13) 2.3 Driven particle Consider a uniform ow past a particle in an unbounded three-dimensional (3D) uid. The solution in this case can be obtained by superimposing the solutions of uniform ow, stokeslet and 3D potential source dipole, where the orientations of the induced stokeslet and potential dipole are pointing opposite to the uniform ow. The ow of a stokeslet decays as 1=r, wherer is the distance from the particle, 22 (b) Active System (a) Driven System Force dipole (Unconfined) Potential dipole (Confined) Stokeslet (Unconfined) Potential dipole (Confined) Figure 2.2: Unconned and conned microswimmers have distinct hydrodynamic signature in far eld in the active and driven system. (a) In the active system, the far-eld in unconned and conned ow induced by swimmer's motion is given by force dipole and potential dipole respectively, (b) whereas in the driven system, the far-eld in unconned and conned ow induced by uniform background ow is given by stokeslet and potential dipole respectively. whereas the ow of a 3D potential dipole decays as 1=r 3 . Therefore, in leading order, the far-eld ow induced by a uniform ow on a particle is given by a stokeslet. However, this far-eld hydrodynamic signature changes drastically in a Hele- Shaw ow. In the conned 2D uid, the leading order ow singularity induced on the particle due to a uniform ow changes from a stokeslet to a 2D potential dipole [10, 73], see Fig. 2.2(a). The scaling of a potential dipole also changes in the conned limit. This change in scaling can be understood by recalling that a potential dipole can be obtained from the modied continuity equation and the velocity potential [6, 13, 33] r u =r(r); u =r; (2.14) 23 where r is the spatial vector and provides the strength and orientation of the potential dipole. The loss of one dimension in the spatial derivative of (2.14) from 3D unbounded uid to the 2D conned system causes the change in scaling of the potential dipole from 1=r 3 decay to 1=r 2 decay. The 2D potential dipole solution is given by u(r) = 1 2jrj 2 (2^ r^ r I); (2.15) where I is the rank-two identity tensor and ^ r = r=jrj. For the case of a particle driven by a uniform ow in the ^ x direction, =a 2 V ^ x, where a is the eective radius of the particle and V is the strength of the uniform ow. Now, we consider a system of N conned driven particles, the long-range inter- actions among the particles can be obtained from the superposition of their induced dipoles [6, 33]. By vitrue of (2.13) and (2.15), for the n-th driven particle, where n = 1;:::;N, its equation of motion is given by _ r n =V ^ x 2 N X j6=n a 2 V 2jr n r j j 2 " 2(r n r j )(r n r j ) jr n r j j 2 I # ^ x: (2.16) 2.4 Swimming particle The far-eld ow induced by the self-propelled motion of a microswimmer can be approximated with multipole expansion [20, 116], that is, using a combination of singularity solutions of the governing equations for the uid ow to approximate the ow induced by the swimmer. The induced ow singularities depend on the 24 microscopic details of the swimmer's propulsion mechanism and the transport of the incompressible uid due to the nite size of the swimmer. In unconned three-dimensional (3D) uid, the far-eld ow induced by a swim- mer is usually approximated by a force dipole (stresslet), arising from the stress distribution on the swimmer's surface that creates the propulsion of the swimmer. The swimmer can be categorized into pusher or puller [102, 103]. A pusher corre- sponds to swimmer that exerts a propulsive force near the rear of its body, whereas a puller, exerts a thrust at its front. The ow direction of the force dipole induced by pusher and puller is exactly opposite to each other, and this dipolar ow decays as 1=r 2 , where r is the distance from the swimmer. In addition to the force dipole, another ow singularity arises from the incompressibility, which is given by the po- tential dipole. However, the potential dipole is subdominant in comparison with the force dipole, and decays as 1=r 3 , hence not contributing to the far-eld ow. Consider the situation when a microswimmer is strongly conned in a Hele-Shaw cell. It turns out that, the scalings of the ow singularities induced by the swim- mer change drastically under the connement, leading to distinct hydrodynamic signature compared to the ow in the unconned case [13, 33]. The ow induced by the force dipole is screened by the conning walls via the shear ow in the thin lubrication lms between the swimmer and the walls, thereby decaying rapidly as 1=r 3 instead of 1=r 2 . This is attributed to the fact that a stokeslet (force monopole) changes its scaling to 1=r 2 in the conned limit [10, 73]. It follows that the force dipole, which can be constructed from the spatial derivative of the stokeslet, decays 25 as 1=r 3 . In contrast to the increase in the order of decay in the ow velocity of force dipole, the potential dipole that arises from incompressibility of the uid decays slowly as 1=r 2 in the conned 2D system, as opposed to the fast 1=r 3 decay in the unconned case. Therefore, in the case of conned microswimmer, the force dipole is a subdominant contribution compared to the potential dipole, and the far-eld ow induced by a swimmer is given by a potential dipole, as opposed to the force dipole in unconned case, see Fig. 2.2(b). The potential dipole induced by a conned microswimmer has the same form as (2.15), with =a 2 U ^ p, where a is the eective radius of the swimmer, U is the self-propelled speed, and ^ p is the orientation of the swimmer. Consider a system ofN conned microswimmer and take into account their long-range hydrodynamic interactions. For the n-th swimmer, where n = 1;:::;N, its translational equation of motion is given by _ r n =U ^ p n + N X j6=n a 2 U 2jr n r j j 2 " 2(r n r j )(r n r j ) jr n r j j 2 I # ^ p j : (2.17) The orientation equations for the swimmers depend on the details of their geometries and will be discussed in detail in Chapter 3. For both driven and swimming particle, the leading order ow singularity is given by potential dipole. However, in the driven system, the potential dipole is always pointing opposite to the background ow, irrespective to the orientation of the particle. To conclude, microswimmers in Hele-Shaw geometries have a distinct hydrody- namic signature in the sense that the far-eld ow is that of a 2D potential source 26 dipole as opposed to the 3D force dipole in the unconned case. Thus, the usual categorization of unbounded swimmers into pushers and pullers becomes irrelevant. The dipolar far eld is independent of the transport mechanism, be it driven parti- cles or self-propelled swimmers, pushers or pullers. It is rooted in the fact that the basic physics is that of a Hele-Shaw ow. 27 Chapter 3 Problem formulation In this chapter, we formulate the potential dipole models for microswimmers in con- ned Hele-Shaw geometries. Our goal is to develop reduced order models that take into account the leading order hydrodynamic interactions among the microswim- mers. These models will be employed in later chapters to study the collective dynamics of large populations of microswimmers in dierent conned micro uidic systems. One main challenge in deriving the dipole models is to obtain the correct orien- tational dynamics of the microswimmers. We show, via simple physical argument, the details of how a swimmer reorient in response to local ow is intrinsically related to the geometry of its body. Here, we consider two categories of idealized minimal models of microswimmers { symmetric and asymmetric dumbbell swimmers, which are composed of two circular discs connected with a frictionless rod. We assume that the dumbbell swimmer self-propels at a constant speed, without considering the details of its propelling mechanism. All these models dier only by the relative 28 size of the two circular discs and the locations of the circular discs relative to the self-propelled motion. For the symmetric swimmers, we consider two dierent locations for the circular discs of the dumbbell, see Fig. 3.1(a). One choice is to place the discs on the line transverse to the self-propelled motion of the dumbbell, which is used to approx- imate a blu body and we refer to this swimmer as the blu swimmer. Another choice would be to place the discs on the line aligned with the self-propelled motion of the dumbbell, which corresponds to a slender body, and it is called the slen- der swimmer. Interestingly, the orientation equations for these two models dier by only a negative sign. Further, it turns out that the slender swimmer model is equivalent to the conned microswimmer model used in [33], and the orientation equation is consistent with the Jeery's equation for the motion of an ellipsoid in shear ow at zero Reynolds number [55]. We also revisit the asymmetric dumbbell swimmer model proposed by Brotto et al. [13]. The asymmetric dumbbell swimmer is characterized by two unequally sized circular discs aligning with the self-propelled motion of the dumbbell, creat- ing a head-tail asymmetry, see Fig. 3.1(b). The head-tail asymmetry causes a given swimmer to reorient, not only in response to the ow gradient as anticipated by Jef- fery's equation, but also in response to the ow velocity itself. This result is rooted in the fact that the lubrication forces between the swimmer and the solid walls hin- der its advection by the uid, inducing unequal translational motility coecients at the swimmer's head and tail. 29 Bluff swimmer Slender swimmer change in orientation (a) Large tail swimmer change in orientation Large head swimmer (b) z 1 z 2 z 1 z 2 z head z tail z head z tail Figure 3.1: Dipole model for conned microswimmers: (a) symmetric (blu and slender) swimmer; (b) asymmetric (large tail and large head) swimmer. The organization of this chapter is as follows: Section 3.1 focuses on the sym- metric swimmer model and we formulate the equations of motions for a system of swimmers. We then focus on the asymmetric swimmer model in Section 3.2. Fi- nally, we present the ow velocity eld of a potential dipole in dierent conned geometries in Section in 3.3. 3.1 Symmetric Microswimmer Model Since Hele-Shaw ow is potential, it is convenient to express the ow eld in concise complex notation. Here, we derive the equations of motion for the symmetric mi- croswimmer. Consider a microswimmer composed of two connected discs of radius R located atz 1 andz 2 respectively, see Fig. 3.1(a), wherez =x + iy is the complex coordinate (i = p 1). Assume the two discs are connected by a frictionless rod 30 of length `. The equations of motion for the two discs z 1 and z 2 can be written in complex notation as _ z 1 =U o e io + w(z 1 ) +T e io ; _ z 2 =U o e io + w(z 2 )T e io : (3.1) Here,U o is the swimmer's self-propelled velocity, o its orientation angle, andw(z) is the velocity eld of the ambient uid. The bar notation denotes the complex conjugate, z = xiy. The coecients is the translational mobility coecients that arises from the balance of hydrodynamic drag and wall friction acting on the discs, which is a decreasing functions ofR. The variable is the inverse of the drag coecient andT is the unknown constraint force, together they form the Lagrange multipliers that enforces the constraintjz 1 z 2 j =`, which is purely imaginary for the blu swimmer and purely real for slender swimmer. The hydrodynamic center of the swimmer is given byz o = (z 1 +z 2 )=2. Our goal is to rewrite the equations of motion (3.1) in terms of the swimmer's hydrodynamic centerz o and orientation o . Assume the complex conjugate ow velocity is analytic, we use Taylor series expansion to expand the ow velocity at the two discs. For simplicity, the self-interaction among the two discs are neglected, as its eect only amounts to a slight modication in the eective translational mobility coecient of 31 the entire swimmer. For the blu swimmer, the ow velocity on the two discs are expanded as w(z 1 ) = w(z o ) i` 2 e io d w dz zo +::: w(z 2 ) = w(z o ) + i` 2 e io d w dz zo +:::: (3.2) and similarly, for the slender swimmer, we have w(z 1 ) = w(z o ) ` 2 e io d w dz zo +::: w(z 2 ) = w(z o ) + ` 2 e io d w dz zo +:::: (3.3) Substitute (3.2) or (3.3) into (3.1), we obtain the equation governing the transla- tional motion of the swimmer's center _ z o =U o e io + w(z o ); (3.4) which has the same form for both blu swimmer and slender swimmer. The dier- ence between the blu swimmer and the slender swimmer arises in their orientation equations, which are derived from the relations governing the location of the two discs z 1 and z 2 . For the blu swimmer, the relation i`e io =z 2 z 1 holds, whereas for the slender swimmer, the location of the two discs are governed by the relation 32 `e io = z 2 z 1 . Upon dierentiating both sides with respect to time and further simplications, we obtain, for the blu swimmer, _ o =Re ( _ z 2 _ z 1 )e io ` ; (3.5) and for the slender swimmer, _ o = Re ( _ z 2 _ z 1 )ie io ` : (3.6) Here, Re denotes the real part of the expression in bracket. Now, we obtain the orientation equation via substituting (3.1) and (3.2) into (3.5) for the case of blu swimmer, and substituting (3.1) and (3.3) into (3.6) for the case of slender swimmer. It turns out that (3.5) for the blu swimmer and (3.6) for the slender swimmer collapse to the same form, except the dierence in sign for the rotational mobility coecient 1 , _ o = 1 Re dw dz ie 2io ; (3.7) where 1 =. The (+) and () signs correspond to the slender and blu swimmer respectively. The orientation equation for the slender body is consistent with the Jeery's equation for the motion of an ellipsoid in shear ow at zero Reynolds number. Note thatj 1 j in general, but since our dumbbell models do not take into account the thickness of the body, 1 is equal to its value at the limit of zero body thickness, therefore,j 1 j =. 33 Now, we consider the interaction of N microswimmers in an unbounded uid domain. Each swimmer induces a far-eld velocity of a potential dipole, [6], the far- eld ow of a microswimmerj located atz j =x j + iy j and oriented at an arbitrary angle j can be described by the complex velocityw(z) =u x iu y =e i j =(zz j ) 2 , where is the dipole strength. Note that =R 2 U, whereR is the eective radius of the swimmer, and U is the self-propelled speed of the swimmers. A microswimmer n responds to the ow induced by all microswimmers in the uid domain, namely, w(z n ) = N X j 6=n j=1 e i j (z n z j ) 2 : (3.8) By virtue of (3.4) and (3.7), the equations governing dynamics of N interacting swimmers, all having the same self-propelled speed U, is given by _ z n =Ue in + N X j6=n e i j (z n z j ) 2 ; _ n = 1 Re " N X j6=n 2ie 2in e i j (z n z j ) 3 # ; (3.9) We conclude this section by writing the system of equations (3.9) in dimen- sionless form using the swimmers radius R as a length scale and R=U as a time scale. That is, we introduce the dimensionless spatial variable ~ z = z=R and time variable ~ t =tU=R. We then drop the tilde notation assuming all variables are non- dimensional. The parameters U and are now normalized to one, that is, U = 1 34 and = 1. The coecients and 1 are also non-dimensional. The dimensionless equations are therefore given by _ z n = e in + N X j6=n e i j (z n z j ) 2 ; _ n = 1 Re " N X j6=n 2ie 2in e i j (z n z j ) 3 # : (3.10) 3.2 Asymmetric Microswimmer Model In this section, we revisit the head-tail asymmetric microswimmer model proposed in [13]) and reformulate the equations of motion in complex notation. Consider a microswimmer composed of two connected discs of radii R tail and R head located at z tail and z head respectively, see Fig. 3.1(b). Assume the two discs are connected by a frictionless rod of length`. The equations of motion for the swimmer's tail (z tail ) and head (z head ) are given by _ z tail =U o e io + tail w(z tail ) + tail T e io ; _ z head =U o e io + head w(z head ) head T e io : (3.11) Here,U o is the swimmer's self-propelled velocity, o its orientation angle, andw(z) is the velocity eld of the ambient uid. The coecients tail , head are the transla- tional mobility coecients whereas tail , head are the inverse of the drag coecients at the tail and the head, and T is the constraint force that enforces the constraint 35 jz head z tail j =`. The translational mobility coecients tail and head arise from the balance of hydrodynamic drag and wall friction acting on the tail and head, and are decreasing functions of R tail and R head respectively, with values less than 1, see, e.g. [13, 33]. We dene the hydrodynamic center of the swimmer to be z o = ( tail z head + head z tail )=( tail + head ). Similar to the symmetric swimmer, we rewrite the equa- tions of motion (3.11) in terms of the swimmer's hydrodynamic center z o and ori- entation o . Let`R tail ;R head , neglect self-induced ow of the two discs, and use Taylor series expansion to expand the ow velocity at the tail and head w(z tail ) = w(z o )`e io tail tail + head d w dz zo +::: w(z head ) = w(z o ) +`e io head tail + head d w dz zo +:::: (3.12) Substitute (3.12) into (3.11) to get that the equation governing the translational motion of the swimmer's center _ z o =U o e io + w(z o ); (3.13) where = ( head tail + tail head )=( head + tail ). To obtain the equation governing the rotational motion of the swimmer, note that, by denition, `e io = z head z tail , which gives, upon dierentiating both sides with respect to time and further simplications, _ o = Re ( _ z head _ z tail )ie io ` : (3.14) 36 Now substitute (3.11) and (3.12) into (3.14) to get _ o = Re 1 dw dz ie 2io + 2 wie io : (3.15) where dw dz and w are evaluated at z o and the constant parameters 1 and 2 are given by 1 = ( head head + tail tail ) head + tail ; 2 = ( head tail )=`: (3.16) The sign of 2 dictates how the swimmer orients in local ow: it aligns to the local ow when 2 > 0, that is, for large tail swimmers for which head tail > 0 (because head=tail is a decreasing function of R head=tail , [6]), and opposite to the local ow when 2 < 0 , that is, for large head swimmers for which head tail < 0. Note that in contrast to the symmetric swimmer, 1 is always positive, and the model for asymmetric swimmer reduces to the model for symmetric slender swimmer when 2 = 0. By virtue of (3.13) and (3.15), the dynamics of N swimmers, all having the same self-propelled velocity U, can be expressed as _ z n =Ue in +w(z n ); _ n = Re 1 dw dz ie 2in + 2 wie in : (3.17) 37 To close the model, one needs to obtain an expression for the uid velocity eld w(z). For the case of a system ofN microswimmers in an unbounded uid domain, the ow velocity eld induced on then-th swimmer is given by (3.8). Using (3.8), we obtain from (3.17) the equations of motion governing the dynamics ofN interacting swimmers _ z n =Ue in + N X j6=n e i j (z n z j ) 2 ; _ n = Re " N X j6=n 1 2ie 2in e i j (z n z j ) 3 + 2 e in e i j (z n z j ) 2 # ; (3.18) where = R 2 U, with R being the eective radius and U being the self-propelled speed of the swimmers. The system of equations (3.18) are now written in dimensionless form. The length scale and time scale used to non-dimensionalize the system are exactly the same as that used in Section 3.1 for the symmetric swimmer, namely, we use the swimmers radius R as a length scale and R=U as a time scale. We introduce the dimensionless spatial variable ~ z =z=R and time variable ~ t =tU=R. We then drop the tilde notation assuming all variables are non-dimensional. Note that (3.18) have the same form after the non-dimensionalization but the parameters U and are now normalized to one, that is, U = 1 and = 1. The parameter values , 1 and 2 are also non-dimensional. 38 3.3 Potential dipole in conned systems Here, we derive the ow eld of a potential dipole in dierent Hele-Shaw geome- tries, including the ow in an unbounded and doubly-periodic domain, and the ow subject to connement of circular or sidewall boundaries. For a potential dipole located at z = z o in an unbounded uid domain, its complex potential is given by F (z) = e io zz o ; (3.19) Here, and o are the strength and orientation of the dipole. The complex con- jugate velocity eld due to the dipole can then be evaluated through w = dF=dz, which gives w(z) = e io (zz o ) 2 : (3.20) The ow of a potential dipole in unbounded domain has a angular symmetry and R w(z)dxdy = 0. The ow eld and the streamlines for a potential dipole are depicted in Fig. 3.2. We then consider the potential dipole in a doubly-periodic domain. A dipole in a doubly-periodic domain has a doubly innite set of images. Thus, evaluatingw(z) requires the evaluation of conditionally-convergent, doubly-innite sums of terms that decay as 1=jzj 2 . It turns out that, the potential dipole in a doubly periodic 39 −5 0 5 −5 0 5 x y −5 0 −5 0 5 5 y x (a) (b) Figure 3.2: (a) Flow velocity eld and (b) streamline of a potential dipole in an unbounded domain. The ow close to the singularity is covered. domain can be written in a closed-form expression in terms of the Weierstrass elliptic function [57, 130], w(z o ) =}(zz o ;! 1 ;! 2 )e in : (3.21) The Weierstrass elliptic function is given by} (z;! 1 ;! 2 ) = 1 z 2 + P k;l 1 (z kl ) 2 1 2 kl , with kl = 2k! 1 + 2l! 2 , k;l2Zf0g, and ! 1 and ! 2 being the half-periods of the doubly-periodic domain. This function has innite numbers of double pole singular- ities located atz = 0 andz = kl , corresponding to the 1=jzj 2 singularities induced by the potential dipoles. The derivation of (3.21) can be found in Appendix A. The ow eld and the streamlines for the potential dipole in a doubly periodic domain are depicted in Fig. 3.3. 40 x y −15 −10 −5 0 5 10 15 −15 −10 −5 0 5 10 15 (a) (b) −15 −10 −5 0 5 10 15 −15 −10 −5 0 5 10 15 x y Figure 3.3: (a) Flow velocity eld and (b) streamline of a potential dipole and its closest images in a doubly-periodic domain. The ow close to the singularity is covered. Now, we consider the potential dipolar ow with boundary connement. The connement doesn't alter the nature of the dipolar ow, but introduces image dipoles that satisfy the impenetrable boundary conditions. We derive the complex conjugate ow eld of a dipole in a circular domain by applying the Milne-Thomson circle theorem [4]. This theorem states that for any irrotational ow with an associ- ated complex potential functionF (z), analytic in the considered domain, a circular boundary of radius R can be constructed at the origin using the modied complex potential function G(z) =F (z) +F (R 2 =z): (3.22) 41 x y −10 −5 0 5 10 −10 −5 0 5 10 (a) (b) −10 −5 0 5 10 −10 −5 0 5 10 x y Figure 3.4: (a) Flow velocity eld and (b) streamline of a potential dipole in a conned circular domain. The circular boundary is denoted in red color. The ow close to the singularity is covered. The second term in the expression enforces the imaginary part of the modied complex potential functionG(z) to be zero atjzj =R, thus creating an impenetrable circular boundary of radius R. Substitute (3.19) into (3.22) and dierentiate with respective to z, we obtain the complex conjugate ow eld of a dipole in a circular domain [131], w(z) = e io (zz o ) 2 + (R 2 =z 2 o )e i(o) (zR 2 =z o ) 2 : (3.23) The second term of (3.23) corresponds to the contribution of the image dipole due to the circular boundary. Therefore, for a dipole located at z o and orientation of o in a circular boundary of radius R, its image dipole is located at R 2 =z o , with 42 orientation o and dipole strength(R 2 =z 2 o ). The ow eld and the streamlines for a potential dipole in the circular connement are depicted in Fig. 3.4. The potential dipole oset from the center of the circular domain no longer holds its angular symmetry, that is, R w(z)dxdy 6= 0, which can be visualized from the deformation of the ow and streamlines close to the circular boundary due to the impenetrability condition. We proceed to derive the ow eld of a potential dipole in a rectangular channel of widthW . The the sidewall connement introduces impenetrability conditions at z =iW=2, which can be accounted for using the method of images. A dipole at z =z o corresponds to two periodic lattices, of period 2iW , of image dipoles located atz o =x o + iy o andz o =x o + i(Wy o ). By virtue of (3.19), the complex potential for a dipole in a singly periodic domain with period 2iW is given by F periodic =e io 1 z + 1 X n=1 1 z 2inW + 1 z + 2inW ! : (3.24) Recall the series representation cotz = 1 z + P 1 n=1 1 zn + 1 z+n . Equation (3.24) can then be rewritten in a closed-form expression, F periodic = 2W e io coth z 2W : (3.25) Therefore, the complex potential for the image system can be written as F = 2W ( coth h 2W (zz o ) i e io + tanh h 2W (zz o ) i e io ) : (3.26) 43 x y −10 −5 0 5 10 −5 0 5 (a) (b) −10 −5 0 5 10 −5 0 5 x y Figure 3.5: (a) Flow velocity eld and (b) streamline of a potential dipole in a rectangular channel. The red line indicates the sidewall boundaries of the channel. The ow close to the singularity is covered. Dierentiate (3.26) with respect to z gives the complex conjugate ow eld [132], w(z) = 2W 2 ( csch 2 h 2W (zz o ) i e io sech 2 h 2W (zz o ) i e io ) : (3.27) Similar to the case of dipole in circular domain, the potential dipole in rectangular channel has no angular symmetry, evidenced by the deformation of the ow and streamlines close to the lateral boundaries, see Fig. 3.5. 44 Chapter 4 Interaction of a pair of conned microswimmers We examine the detailed dynamics of a pair of symmetric and asymmetric mi- croswimmers introduced in Chapter 3. The symmetric swimmer can be distin- guished into blu or slender swimmer, characterized by two connected circular discs placing at the locations transverse to or along the self-propelled motion of the swimmer, see Fig. 3.1(a). The asymmetric swimmer is characterized by a head- tail asymmetry, which can be distinguished into large tail or large head swimmer, see Fig. 3.1(b). Here, we focus on the interaction of a pair of microswimmers in unbounded or doubly-periodic domain, which will serve to shed light on the physics governing the collective dynamics of populations of microswimmers. We rst consider the interactions of a pair of blu or slender swimmers in an unbounded domain, comparing and contrasting their behavior as a function of initial conditions. We nd that they exhibit regular behavior consisting of either colliding, synchronizing or diverging trajectories depending on their relative initial spacing 45 and relative initial orientation. We then consider the interactions of a pair of large tail or large head swimmers in a doubly-periodic domain and focus on a special class of solutions where the two swimmers move with constant speed and at constant orientation. We nd two families of these \equilibrium-like" solutions: (1) both swimmers swim side by side in a synchronized way; and (2) one swimmer tailgates the other. We analyze their stabilities and nd that they depend on the details of the head-tail asymmetry. Finally, we discuss the results of interaction of swimmer pairs to highlight their insight on the large-scale simulations. 4.1 Synchronization, divergence and collision of symmetric swimmers Consider a complex plane z = x + iy, and i = p 1. We examine the behavior a pair of symmetric swimmers in an unbounded domain, with locations and orienta- tions denoted by (z 1 ; 1 ) and (z 2 ; 2 ). The equations of motions for two symmetric swimmers, in dimensionless form, are given by _ z 1 = e i 1 + e i 2 (z 1 z 2 ) 2 ; _ 1 = 1 Re 2ie 2i 1 e i 2 (z 1 z 2 ) 3 ; _ z 2 = e i 2 + e i 1 (z 2 z 1 ) 2 ; _ 2 = 1 Re 2ie 2i 2 e i 1 (z 2 z 1 ) 3 : (4.1) The behavior of the swimmer pair is examined by numerically integrating these equations with dierent initial conditions. Sample trajectories are depicted in 46 0 30 60 90 120 −10 −5 0 5 10 15 x y 0 20 40 60 0 2 4 6 8 x y 0 5 10 15 20 0 1 2 3 4 5 x y (a) (b) (c) Figure 4.1: Blu swimmer pairs exhibit (a) collision, (b) synchronization and (c) divergence of their paths. Parameter values: = 0:5, 1 =0:5, z 1 (0) = 1 (0) = 2 (0) = 0. (a) z 2 (0) = 5i; (b) z 2 (0) = 1 + 2:5i; and (c) z 2 (0) = 5 + 2:5i. Total integration time: (a) 50; (b) 60; (c) 100. The positions of the dipoles are marked by `' at t = 0 and by `o' at the end of the integration time. Fig. 4.1 for the blu swimmers and Fig. 4.2 for the slender swimmers. The blu swimmer aord three modes of interactions depending on choice of initial conditions. The three modes are collision, synchronization or divergence. Fig. 4.1(a) shows an example of the two swimmers approaching collision without colliding in nite time. Fig. 4.1(b) shows the synchronization mode where the two swimmers trace peri- odic trajectories with anti-phase symmetry. Fig. 4.1(c) shows an example where the swimmers' paths diverge to innity. Interestingly, the slender swimmers exhibit two modes of interaction only: synchronization or divergence as shown in Fig. 4.2. Poincar e sections for the synchronized trajectories of Fig. 4.1(b) and Fig. 4.2(a) are depicted in Figs. 4.3(a) and 4.3(b), respectively. Note that in both cases, the orientations 1 and 2 of the two swimmers are related by a phase shift and thus it suces to show section (z 1 z 2 ; 1 ) since section (z 1 z 2 ; 2 ) is identical. 47 0 10 20 30 −10 −5 0 5 10 15 x y 0 20 40 60 80 100 −2.5 −1.5 −0.5 0.5 1.5 2.5 x y (a) (b) Figure 4.2: Slender swimmer pairs exhibit (a) synchronization and (b) divergence of their paths. Parameter values: = 0:5, 1 = 0:5,z 1 (0) = 1 (0) = 2 (0) = 0. (a) z 2 (0) = 2 + 2i; (b) z 2 (0) = 2i. Total integration time: (a) 90; (b) 30. The positions of the dipoles are marked by `' at t = 0 and by `o' at the end of the integration time. −1.1 −1.05 −1 −0.95 −2.5 −2 −1.5 −1 −0.2 0 0.2 y 1 − y 2 x 1 − x 2 α −2.2 −2.1 −2 −1.9 −1.8 −2 −1 0 1 2 3 −0.5 0 0.5 y 1 − y 2 x 1 − x 2 α (a) (b) Figure 4.3: Poincar e section for (a) the case shown in Fig. 4.1(b) and (b) the case shown in Fig. 4.2(a). Note that the Poincar e sections (x 1 x 2 ;y 1 y 2 ; 1 ) and (x 1 x 2 ;y 1 y 2 ; 2 ) are coincident because 1 and 2 for these trajectories dier by a phase shift only. We proceed to investigate the dependence of the dierent modes of behavior on initial conditions. The results are depicted in Fig. 4.4. We let z 1 (0) = 0 + i0, 1 (0) = 0 and varyz 2 (0) and 2 (0) such that, for each value of 2 (0),z 2 (0) is varied in a square box from -20 to 20 (dimensionless units) using xed small increments 48 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 (a) (b) (c) (d) (e) −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 −10 −5 0 5 10 Figure 4.4: Blu (top row) and slender (bottom row) swimmer: phase diagrams showing the behavior of the blu and slender swimmer pair as function of initial conditions. The horizontal and vertical axis of the diagrams correspond to the real and imaginary part of the initial relative spacing of the swimmers z 2 (0) z 1 (0) respectively. The black regions correspond to initial conditions where the two swimmers diverge, the grey regions correspond to synchronized motion of the two swimmers, and the white regions (except the excluded white disc around the origin) correspond to colliding trajectories. Note that x and y in the plots are in units of p . x + iy. The interaction mode (collision, synchronization, or divergence) is then recorded as a function of the initial spacing of the swimmers. Regions of distinct be- havior are depicted in dierent colors: white for collision, grey for synchronization, and black for divergence. To better resolve the boundaries between these regions, a shooting method is used to rene the initial positionz 2 (0) at the boundary. Note that the small white circle with radius of 1 around the origin is excluded from this study because it corresponds to a very small distance between the two swimmers less than one body length of a swimmer. Figs. 4.4(a) { (e) correspond to the phase diagrams of blu swimmer pair. Clearly, the synchronization and collision regions 49 decrease as the initial relative orientation 2 (0) 1 (0) between the two swimmers increases. To explain this, one can argue intuitively that two swimmers facing oppo- site directions would be more likely to diverge than synchronize or collide together. Another important observation is that, for all values of 2 (0), the diagram (except the collision region) possesses a re ection symmetry about two axes (e i(0) ; ie i(0) ), where is dened as = ( 1 + 2 )=2. Figs. 4.4(f) { (j) correspond to the phase diagrams of slender swimmer pair. No collision behavior is observed. The space of initial conditions is split into synchro- nization and divergence regions. As in the case of the blu swimmers, the syn- chronization region decreases in size as the initial relative orientation 2 (0) 1 (0) between the two swimmers increases. Here again, the diagram possesses a re ec- tion symmetry about the two axes (e i(0) ; ie i(0) ). This re ection symmetry is not a coincidence, but indeed a fundamental property of the two swimmers system. More generally, the behavior of these swimmer pairs, when expressed in terms of a frame moving with = ( 1 + 2 )=2, can be described by the relative spacing and relative orientation of the swimmer pair, as we show next. To this end, let = ( 1 2 )=2 and consider the transformation of coordinates from the inertial complex plane to a plane co-rotating with such that z () = () e i ; e i 1 = e i e i ; e i 2 = e i e i : (4.2) 50 Substitute into (4.1) and simplify, using = 1 2 and = 1 + 2 , to get _ = i _ + 2 cos 1 + 2 ; _ =2 1 Im 1 3 cos; (4.3) and _ = i _ 2i sin 1 + 2 ; _ = 2 1 Re 1 3 sin: (4.4) One can readily verify upon substituting the expression for _ from (4.3) into (4.4) that (4.4) gives an independent system of equations for (;) that can be solved separately. Borrowing terminology from geometric mechanics, the \shape" vari- ables (;) determine the relative spacing and relative orientation of the swimmer pair whereas the variables (;) determine its location and orientation relative to the inertial frame and can be \reconstructed" posteriori from (4.3) once (;) are computed using (4.4). This reasoning is also applicable in the general case ofN swimmers. There, the reduced 3(N 1) shape equations are obtained by a transformation of coordinates from the inertial frame to a frame co-rotating with = P N n=1 n =N and where = P N n=1 z n =N. The equations for these variables should decouple from the rest of the system. An interesting question that remains to be answered is whether the interactions of N swimmers conserve any integrals of motion and/or possess a Hamiltonian structure. Certainly factoring out the equations governing the position and orientation of the system is possible because the system has translational 51 and rotational symmetries. However, the conserved quantities associated with these symmetries are not easily identiable and remain the topic of a future study. 4.2 Pursuit and synchronization of asymmetric swimmers We consider two asymmetric swimmers in a doubly-periodic domain. The locations and orientations of the swimmers are denoted by (z 1 ; 1 ) and (z 2 ; 2 ). The equation of motions of such swimmers, in dimensionless form, are given by _ z 1 = e i 1 +}(z 1 z 2 )e i 2 ; _ z 2 = e i 2 +}(z 2 z 1 )e i 1 ; _ 1 = Re 1 } 0 (z 1 z 2 )ie i(2 1 + 2 ) + 2 }(z 1 z 2 )ie i( 1 + 2 ) ; _ 2 = Re 1 } 0 (z 2 z 1 )ie i(2 2 + 1 ) + 2 }(z 2 z 1 )ie i( 2 + 1 ) : (4.5) Here,}(z;! 1 ;! 2 ) is the Weierstrass elliptic function, which is given by} (z;! 1 ;! 2 ) = 1 z 2 + P k;l 1 (z kl ) 2 1 2 kl , with kl = 2k! 1 + 2l! 2 , k;l2 Zf0g, and ! 1 and ! 2 being the half-periods of the doubly-periodic domain. This function has innite numbers of double pole singularities located atz = 0 andz = kl , corresponding to 52 0 π/4 π/2 3π/4 π 0 π/4 π/2 3π/4 π 5π/4 3π/2 α θ stable unstable unstable stable pursuit synchronization large tail large head ν > 0 ν < 0 (a) (b) Stability Analysis α1 α2 θ c 2 2 Figure 4.5: (a) Two modes of solutions are obtained from (4.7): a pursuit mode (blue) where the swimmers trail one another and a synchronization mode (red) where they swim side by side. The shown solid curves correspond to! 1 =i! 2 = 5 and the dashed lines to ! 1 =i! 2 = 10. The separation distance c is set to c = 4. (b) Summary of the stability analysis for these two modes. the 1=jzj 2 singularities induced by the potential dipoles. By introducing this func- tion, we take into account the interaction among the swimmers and their periodic images. We focus on the dynamic response of the swimmers when 1 = 0, that is, when their alignment with the ow gradient is negligible. In this case, the orientation dynamics is dominated by alignment with the ow due to head-tail hydrodynamic asymmetry. Since the alignment to local ow gradient is a higher order eect compared to the alignment to local ow in (4.5), the approximation of 1 = 0 is valid when the swimmer is suciently asymmetric, that is, the value of 2 is nite. We look for special solutions where the two swimmers move at the same velocity and orientation for all time. That is, we look for solutions where _ z 1 = _ z 2 = constant 53 and _ 1 = _ 2 = 0. To obtain the initial conditions that lead to this behavior, it is convenient to rewrite the equations of motion (4.5) in terms of the reduced coordinate z 1 z 2 which we set to z 1 z 2 = =ce i (see inset of Fig. 4.5(a)). To this end, one gets _ = e i 1 e i 2 +}()(e i 2 e i 1 ); _ 1 = _ 2 = 2 Re[ie i( 1 + 2 ) }()]: (4.6) The translation equation for _ is identically zero when 1 = 2 = . Whereas to guarantee _ 1 = _ 2 = 0, one must satisfy the condition 8 > > > > > > > > > > < > > > > > > > > > > : Re[}(ce i )] = 0 therefore = 4 ; 3 4 and (;) = 8 > < > : (;) (; +=2) or Im[}(ce i )] Re[}(ce i )] = tan 2; Re[}(ce i )]6= 0; (4.7) A total of ten strict relative equilibria of the two swimmers are depicted schemat- ically in Fig. 4.5(a). Namely, the ve solutions given by = 0; 4 ; 2 ; 3 4 ; and = correspond to the two dipoles moving parallel to each other in a \pursuit" mode (blue arrows), whereas the ve solutions given by = 0; 4 ; 2 ; 3 4 ; and = +=2 correspond to the two dipoles moving in tandem in a \synchronized" mode (red arrows). 54 x -10 -5 0 5 10 y -10 -5 0 5 10 x -10 -5 0 5 10 y -10 -5 0 5 10 (a) (b) Figure 4.6: Aperiodic behavior of two dipoles in doubly-periodic domain. The parameter values are = =3, z 1 (0) =2 exp(i), z 2 (0) = 2 exp(i), while is obtained by numerically solving the rst condition in (4.7). The positions of the dipoles are marked by `' att = 0 and by `o' at the end of the integration. As time progresses, the two trajectories densely ll the whole domain. The second set of solutions, that is, the values of (;) for which the conditions of solutions are given by Im[}(ce i )]=Re[}(ce i )] = tan 2 and Re[}(ce i )]6= 0, are not analytically available and need to be computed numerically. Fig. 4.5(a) shows the values of (;) that satisfy these conditions { clearly, two branches of so- lutions are obtained. These solutions depend implicitly on the domain size (! 1 ;! 2 ) and on c, the separation distance between the two swimmers. In other words, for a choice of domain size and separation distance c, (;) are computed accordingly such that the two dipoles move at the same constant velocity and orientation for all time. The two branches shown in Fig. 4.5(a) correspond to two modes of behavior: a pursuit mode where one swimmer trails the other, and a synchronization mode 55 where the two swimmers move side by side. These solutions, while they correspond to the dipoles moving at constant velocity and orientation, can exhibit two dis- tinct types of dynamical behavior due to the doubly-periodic nature of the domain, namely, they could lead to aperiodic and periodic motion of the dipoles. Aperi- odic motion refers to the case where the paths of the dipoles densely ll the whole domain, as shown in Fig. 4.6. This seems to be the generic behavior for arbitrary initial conditions. Periodic behavior refers to trajectories that satisfy the condition z 1 (T ) = z 1 (0) + 2p! 1 + 2q! 2 ; z 2 (T ) = z 2 (0) + 2p! 1 + 2q! 2 ; (4.8) where p and q are integers and T is the period of the motion. This amount to the additional condition 1 (0) = 2 (0) = tan 1 ( q p ): (4.9) The ratio of q=p indicates the ratio of the number of times the dipole crosses the y andx axes in one periodT . Fig. 4.7 depicts the periodic behavior of two dipoles in pursuit and synchronization modes for q=p = 3. We analyze the linear stability of the pursuit and synchronization modes by considering small perturbations = x + i y , 1 and 2 about = ce i ( x = c cos, and y = c sin) and 1 = 2 = , with (;) satisfying (4.7). We 56 −10 −5 0 5 10 −10 −5 0 5 10 x y −10 −5 0 5 10 −10 −5 0 5 10 x y (a) (b) Figure 4.7: (a) Pursuit and (b) synchronization in two dipoles undergoing periodic motion. The parameter values are = tan 1 (3), z 1 (0) = 2 exp(i), z 2 (0) = 2 exp(i), while is obtained by numerically solving the rst condition in (4.7). The positions of the dipoles are marked by `' att = 0 and by `o' at the end of the integration time t = 80. linearize equations (4.6) accordingly. The linearized equations can be written in matrix form as follows: d dt 0 B B B B B B B B @ x y 1 2 1 C C C C C C C C A =M 0 B B B B B B B B @ x y 1 2 1 C C C C C C C C A ; (4.10) 57 −10 −5 0 5 10 −10 −5 0 5 10 x y −10 −5 0 5 10 −10 −5 0 5 10 x y (a) (b) Figure 4.8: Pursuit and synchronization modes as attracting modes. Two dipoles hone in on quasi periodic trajectories where: (a) one dipole is in pursuit of the other for 2 = 0:5. (b) the two dipoles synchronize and move along parallel trajectories for 2 =0:5. The initial conditions are z 1 (0) = 1 (0) = 2 (0) = 0 while z 2 (0) = 1:5 + 1i in (a) and z 2 (0) = 2 + 1i in (b). where the Jacobian matrix M is given by M = 0 B B B B B B B B @ 0 0 sinRe[ie i }()] sin+Re[ie i }()] 0 0 cosIm[ie i }()] cos+Im[ie i }()] 2 Re[ie i2 } 0 ()] 2 Re[e i2 } 0 ()] 2 Re[e i2 }()] 2 Re[e i2 }()] 2 Re[ie i2 } 0 ()] 2 Re[e i2 } 0 ()] 2 Re[e i2 }()] 2 Re[e i2 }()] 1 C C C C C C C C A : (4.11) We compute the eigenvalues numerically and nd that, for large tail swimmers 2 > 0, the pursuit mode is stable, whereas for large head swimmers 2 < 0, the synchronization mode is stable. Our ndings are summarized in Fig. 4.5(b). 58 We test our results numerically by integrating the nonlinear equations (4.5) for arbitrary choices of initial conditions. Interestingly, the pursuit and synchronization modes seem to be globally attracting modes in the case of large tail and large head swimmers, respectively. Fig. 4.8(a) shows a depiction of two large-tail swimmers honing in on quasi-periodic pursuit trajectories, while (b) depicts two large-head swimmers synchronizing their motion in nite time to swim side by side. We conclude this section by noting that in the limit of innite domain, the solutions (4.7) of the doubly-periodic system (4.5) converge to the relative equilibria of the unbounded system. In the unbounded system, the relative equilibria can be obtained either by symmetry arguments or by analytical manipulation of the equations of motion. Namely, one has two families of relative equilibria = and = +=2 which correspond to pursuit and synchronization trajectories, respectively. The convergence of the solutions in (4.7) to these solutions is relatively fast, as indicated in Fig. 4.5(a). As! 1 ,! 2 !1, the Jacobian matrixM converges to M 1 = 0 B B B B B B B B B B @ 0 0 sin+ c 2 sin(2) sin c 2 sin(2) 0 0 cos c 2 cos(2) cos+ c 2 cos(2) 2 2 c 3 sin(23) 2 2 c 3 cos(23) 2 c 2 cos(22) 2 c 2 cos(22) 2 2 c 3 sin(23) 2 2 c 3 cos(23) 2 c 2 cos(22) 2 c 2 cos(22) 1 C C C C C C C C C C A : (4.12) The corresponding eigenvalues are [0, 0, 0, 2 2 =c 2 ]. The eigenvalue2 2 =c 2 corresponds to the pursuit mode where =, whereas +2 2 =c 2 corresponds to the 59 synchronization mode. This means that, for large tail swimmers with 2 > 0, the pursuit mode is linearly stable and the synchronization mode is unstable, whereas for large head swimmers when 2 > 0, the opposite is true, thus conrming the results obtained above for nite-sized doubly-periodic domains. 4.3 Discussion In this chapter, we considered the hydrodynamic dipole models governing the in- teraction of symmetric and asymmetric microswimmers in Hele-Shaw connement, and treated in details the dynamics of two interacting swimmers. For symmetric swimmers, we observed three distinct modes of behaviors, including collision, di- vergence and synchronization in the case of blu swimmer, whereas for the slender swimmer, only the divergence and synchronization modes are observed. For asym- metric swimmers, we found two special solutions that correspond to pursuit and synchronization of the two swimmers. The pursuit mode is stable and attracting for large tail swimmers while the synchronization mode is stable and attracting for large head swimmers. By attracting, we mean that, starting from arbitrary initial conditions, large tail swimmers tend to tailgate each other while large head swimmers tend to synchronize their motion in nite time to swim side by side. Note that the results of interaction between asymmetric swimmer pair are par- ticularly interesting in light of the collective behavior exhibited by populations of such swimmers. The stable pursuit mode for large tail swimmer pair suggests that large tail swimmers tends to tail-gate each other, and it turns out that this property 60 persists for the swirling motions for a population of large-tail swimmers, which will be shown in Chapter 5. It is also worth to highlight that the stable synchronization modes for large head swimmers may be responsible for the emergence of polarized density waves by populations of large head swimmers observed in [69], in which the simulations are performed with a desingularized potential dipole model. How- ever, this polarization mode is not observed in [130] and Chapter 5 which use a local repulsion potential for collision avoidance and accurately account for hydro- dynamic interactions and the doubly innite image system. While the dierence in numerical implementation may play a role, the system size may be the main reason why polarized waves are not observed in [130]. Brotto et al. (2013) predicted this behavior for a continuous kinetic-like model; therefore, it is not surprising that it is not reproduced by a fully nonlinear model with only a few hundred swimmers. Irrespective of the reason, the results reported in this study suggest that the global polarization mode in large-head swimmers is not robust to system disturbances, whereas the pursuit mode in large-tail swimmers is. 61 Chapter 5 Flagella-induced transition in the collective behavior of conned microswimmers 5.1 Introduction Collective motion emerges in a wide range of natural systems, from sh schools [25, 26] to bacterial colonies [21, 28, 36, 114, 115, 148], and is believed to play important roles in the functioning and survival of the group. For example, many species of bacteria cyclically transition from a free-swimming planktonic state into a sessile biolm state in response to environmental conditions [37, 82]. It is well known that the free-swimming state is characterized by agella-driven motility, which is suppressed in the biolm state. However, the eect of intermediate states { where variability in the environmental conditions and/or energy supply to the agellar 62 propulsive system may alter the level of agellar activity { on the emergent bacterial behavior remains largely unexplored. Geometric connement is a common feature of the natural environment of var- ious bacteria species. It is also characteristic of several experimental set-ups on biological and articial microswimmers [12, 28, 126, 148]. Conned microswimmers have a distinct hydrodynamic signature in the sense that the far-eld ow is that of a 2D potential source dipole. Further, due to friction with the nearby walls, conned swimmers with geometric polarity (large head or large tail) reorient in response to both the local ow eld and its gradient [13], and their collective behavior exhibit instabilities that are qualitatively distinct from the ones observed in unbounded swimmers [13, 69]. In this chapter, we modify the conned asymmetric (head-tail) swimmer model introduced in Chapter 3 and [13] to take into account the eect of agellar activity and steric interaction. Using this modied model, we show that decreasing agellar activity could lead microswimmers, via hydrodynamic interactions only, to transi- tion from a turbulent-like swimming state, akin to the one observed in numerous experiments, to clustering and aggregation. It is also worth noting that asymmetric self-propelled particles interacting sterically, not hydrodynamically, were recently shown to exhibit rich phenomena including aggregate formation [140]. 63 5.2 Modeling agellar activity To formulate a particle swimmer model that takes into account the eect of agellar activity, we rst establish that the intensity of the dipolar far-eld induced by a beating agellum depends on the level of agellar activity: vigorously-beating agella induce stronger dipolar elds than weakly-beating ones. We consider a gap-averaged two-dimensional beating agellum in potential ow, governed by the Laplace's equation for ow velocity potential. The gap-averaged agellar motion is prescribed as y(x;t) = A cos(kxt), with x2 [1; 1], and is assumed to induce a constant swimming velocity U in thex-direction. Here, all parameters are dimensionless with the characteristic length and time scales being set by the agellum's half-length and beating frequency, respectively. For simplicity of calculation, we choose the swimming frame that is moving with the beating agel- lum at velocityU as the frame of reference [68], see Fig. 5.1. For small deformations, the transverse ow velocity on the beating agellum is given by @y=@t +U@y=@x, which provides the boundary condition for the problem. The potential ow pertur- bation induced by this beating motion is computed numerically using panel method. We compute the generated ow for various values of the beating amplitude A and wavelength k while normalizing the swimming velocity U to 1. We average the induced ow eld over one beating cycle. The time averaged ow eld is then approximated, using a standard tting method based on minimization of theL 2 -norm, by the dipolar eld of a circular disk with eective radiusR tail moving at the same swimming velocity U, see Fig. 5.2. That is to say, we minimize the 64 −U u = 0 ∞ 1 / k − U Laboratory frame u = U ∞ 1 / k Swimming frame Figure 5.1: Kinematic for the planar waving motion of a beating agellum in the laboratory frame and the swimming frame. Fundamental singularities of potential ow, illustrated by the black dots, are distributed along the agellum to calculate the ow induced by its beating motion. function R (jw tail j 2 jw disc j 2 )dxdy, where w tail is the averaged ow generated by the beating agellum andw disc is the dipolar ow induced by a circular disc of radius R tail . Note that the dipolar eld induced by a circular disk located at z o =x o +iy o (i = p 1) in the complex z-plane and oriented at an arbitrary angle o to the x-axis can be described by the complex velocity w(z) =u x iu y =e io =(zz o ) 2 , where = UR 2 tail is the dipole strength. To this end, we set 0 = and vary z 0 and R tail to minimize the above L 2 -norm. Fig. 5.2(c) shows that, as A and k increase, R tail increases accordingly. Consequently, the dipole strength increases with increasing agellar activity. We now model the agellar far-eld ow by that of a circular disk of eective radius R tail and we assume that the hydrodynamic-coupling between the head and the agellum is weak. This leads to a head-tail dumbbell swimmer model, where the value ofR tail is interpreted as a measure of the agellar activity, see Fig. 5.2(d). We can now match the agella-based head-tail dumbbell swimmer model with the asymmetric dumbbell swimmer model presented in previous chapters: a swimmer with vigorously beating agella is analogous to a large tail swimmer and tends to 65 − 8 −4 0 4 8 −8 −4 0 4 8 − 8 −4 0 4 8 −8 −4 0 4 8 weakly coupled (a) (c) (d) R tail Head Flagellum 1 2 3 4 5 0 0.5 1 1.5 k(π) R A = 1 A = 0.5 A = 0.1 tail R head μ head μ tail (b) t = π/2 t = π t = 3π/2 Figure 5.2: (a) Time-averaged ow eld created by a beating agellum withA = 0:5, k = and U = 1. (Inset) Snapshots of the agellum-induced ow eld at dierent times. (b) Dipolar eld created by a circular disk of eective radius R tail tted to the average ow eld in (a). (c) Change inR tail withA andk of the traveling wave via the agellum. (d) Reduction of a agellated swimmer to a dumbbell swimmer. 66 align to the local ow, whereas a swimmer with weakly beating agella is analogous to a large head swimmer and tends to align to the opposite of the local ow. 5.3 Hydrodynamic model for conned microswim- mers in doubly-periodic domain We propose an idealized physical model that based on the head-tail asymmetric dumbbell model introduced in Chapter 3 and [13] to investigate the eects of ag- ellar activity on the hydrodynamic interactions among a population of microswim- mers. The model considers microswimmers that are strongly conned in a thin lm of Newtonian uid, with the dimension of the swimmers being comparable to the thickness of the uid lm, see Figure 5.3. In particular, we consider the mi- croswimmers in a doubly-periodic domain. The dynamics of a population ofN such swimmers can be expressed in concise complex notation _ z n =Ue in +w(z n ) +V n ; _ n = Re 1 dw dz ie 2in + 2 wie in : (5.1) Here,z n and n denote the position and orientation of each swimmer (n = 1;:::;N), Re denotes the real part of the expression in bracket, whereas, 1 and 2 are non- dimensional translational and rotational mobility coecients whose values depend 67 L L h h y x α Figure 5.3: Strongly conned microswimmers in a Hele-Shaw cell create poten- tial dipolar far-eld ows. Microswimmers with uniformly randomly oriented and spatially homogenous distribution in a doubly periodic domain are considered. on the translational mobility coecients head and tail at the swimmer's head and tail. A swimmer can have dierent mobility coecients at its head, denoted by head and at its tail, denoted by tail , with tail being a decreasing function ofR tail . Given that head is constant, the value of head tail can thus change sign from positive to negative by decreasing the level of agellar activity, and vice versa. For example, according to an experimental study on Rhizobium lupini [109], the average head size corresponds to R head 0:4 m whereas depending on environmental conditions (viscosity, temperature, pH value), the agella-beating amplitude is reported to vary between 0 and 0:4 m. Now, assuming an average wavenumber of 5, this variation in amplitude corresponds, in our model, to a variation in R tail from 0 to 0:62m, thus causing a change in sign ofR head R tail and hence head tail . In the 68 model, this signed dierence dictates the signed value of 2 ( head tail ), and therefore how a swimmer reorients in response to the local ow. Vigorous agellar activity for which head tail > 0 (i.e., 2 > 0) causes swimmers to reorient in the direction of the local ow whereas swimmers with weakly-beating agella for which head tail < 0 (i.e., 2 < 0) reorient in the opposite direction to the local ow. Swimmers also reorient in response to the ow gradient as indicated by the 1 term in (6.1), consistently with the classical Jeery's orbit [55]. Note that 1 is proportional to ( head + tail )=2. To close the model in (5.1), we need to evaluate the velocity eld w(z) induced byN potential dipoles in a doubly-periodic domain, which involves the evaluation of conditionally-convergent, doubly-innite sums of terms that decay as 1=jzj 2 . The velocity eld induced by N potential dipoles located at z n with orientation n , n = 1;:::;N, in a doubly-periodic domain can be written in terms of the Weierstrass elliptic function [57, 130] w = N X n=1 n }(zz n ;! 1 ;! 2 )e in : (5.2) Here, n = R 2 U is the strength of the potential dipole associated with the n th swimmer, withR being the eective radius for the whole swimmer. The Weierstrass elliptic function (z) is given by } (z;! 1 ;! 2 ) = 1 z 2 + P k;l 1 (z kl ) 2 1 2 kl , with kl = 2k! 1 +2l! 2 ,k;l2Zf0g,and! 1 and! 2 being the half-periods of the doubly- periodic domain. This function has innite numbers of double pole singularities 69 located at positions of z = 0 and z = kl , corresponding to the 1=jzj 2 singularities induced by the potential dipoles. In addition to the hydrodynamic coupling, we account for steric interactions in (5.1) using a collision avoidance mechanism. That is, we introduce a near-eld repulsive velocity V n based on the repulsive part of the Leonard-Jones potential V n =f 0 N X j6=n 2R jz n z j j 13 z n z j jz n z j j : (5.3) f 0 is a scaling parameter that characterizes the repulsion strength. These near eld interactions decay rapidly outside a small excluded area centered around z n . Their rapid decay ensures that the order of the far-eld hydrodynamic interactions is preserved. 5.4 Emergence behaviors: chaotic swirling, global orientational order, and stable aggregation We focus on the evolution of populations of micro-swimmers that are initially ran- domly oriented but spatially homogenous. We normalize U and n to 1 and set f 0 = 1, = 0:5, N = 100. We vary 1 , 2 and the area fraction A = NA=A 0 , where A is the area of the micro-swimmer, and A 0 =L 2 is the size of the doubly- periodic square domain. We perform Monte-Carlo type simulations in the sense that, for each set of parameters ( 1 , 2 , A ), we run multiple trials corresponding 70 0 0.5 1 1.5 2 0 0.05 0.1 0.15 0.2 velocity Probability −25 −15 −5 5 15 25 −25 −15 −5 5 15 25 0 0.5 1 1.5 0 0.2 0.4 0.6 0.8 velocity 25 −15 −10 −5 0 5 10 15 15 10 5 0 −5 −10 −15 −25 −15 −5 5 15 25 −25 −15 −5 5 15 25 0 0.5 1 1.5 2 0 0.1 0.2 0.3 0.4 velocity 2 = -1 ν 1 = 0.3 Φ A = 0.1 2 = 3 ν 1 = 0.3 Φ A = 0.3 2 = 0 ν 1 = 0.3 Φ A = 0.1 ν ν ν (a) (b) (c) Figure 5.4: (a) Swirling motion in vigorously-beating agellated swimmers; (b) orientational order; (c) aggregation in weakly-beating agellated swimmers. Prob- ability distributions of velocity of the swimmers are depicted in the lower row. to dierent samples of initial conditions. Each sample corresponds to a set of initial orientations n (0),n = 1;:::;N, taken from a uniform probability distribution func- tion ranging from 0 to 2. To ensure spatial homogeneity, the swimmers' centers z n (0) are placed on a regular meshz kl =kx+ily of size x = y =L= p N, but with a superimposed perturbation that randomizes z n (0) around z kl . This pertur- bation is sampled from a uniform probability distribution centered at zero ranging fromL=(2 p N) to L=(2 p N). We observe the emergence of three distinct types of global structures: swirling behavior for vigorously-beating agella, orientational order for a narrow range of parameter values around 2 = 0, and aggregation or clustering for weakly-beating 71 (a) (b) Figure 5.5: (a) Swimmers with vigorously-beating agella orient with local ow and thus tend to chase each other; (b) Swimmers with weakly-beating agella orient in the opposite direction to local ow and thus tend to aggregate together. agella. Representative simulations are shown in Fig. 5.4. The swirling-like motion is characterized by a velocity distribution function with a mean value higher than the speed of the individual swimmer (Fig. 5.4(a)). This collective behavior where vortex-like structures emerge, break, and rearrange elsewhere is reminiscent to the bacterial turbulence observed in numerous experiments [21, 28, 34, 36, 81, 114, 115]. The associated increase in swimming speed was also observed experimentally. In the context of the dipole model, the increase in speed can be explained as follows. Swimmers with vigorously-beating agella tend to \tail-gate" each other as a result of them aligning with the local ow eld. When a swimmer is close to another swimmer, it will orient towards and travel along the streamlines of the potential dipole created by the nearby swimmer as depicted schematically in Fig. 5.5(a). As the swimmers align and form a chain-like structure, they create a ow eld that helps their forward motion, thus increasing the swimmers' velocities, as evidenced from the probability distribution function in Fig. 5.4(a). 72 As the agellar activity decreases (by decreasing the value of 2 ), a transitional behavior that does not exhibit a clear pattern is observed before a global orien- tational order develops, see Fig. 5.4(b). The development of orientational order happens at a much longer time scale than the swirling dynamics and is a result of the dipoles reorienting with the local velocity gradient (non-zero 1 ), which is ignored in [13, 69]. In this case, the velocity distribution function is Gaussian cen- tered at U = 1 (Fig. 5.4(b)), implying that this collective mode has no advantage at the population level in terms of increased swimming speed. When we further decrease the agellar activity, the swimmers begin to aggregate and cluster, see Fig. 5.4(c). The clustering behavior of swimmers with weak agellar activity can be explained by recalling that such swimmers reorient in the opposite direction of the local ow eld. Therefore, a swimmer tends to travel in the opposite direction to the streamlines created by a nearby swimmer, see Fig. 5.5(b), which leads to aggregation. The collective aggregation of many swimmers takes place at a very short time scale and leads to the formation of clusters. Some clusters are unstable and break readily (blue circle in Fig. 5.4(c)). However, long-lived stable clusters also form (red circle in Fig. 5.4(c)). The stable cluster depicted in Fig. 5.4(c) slows down considerably as more swimmers join the cluster. The velocity distribution function has a strong peak around zero velocity because most swimmers are attached to the stable cluster, with the remaining swimmers moving at unit speed. 73 5.5 Statistics of global modes We use a number of statistical measures to assess the observed global structures. In addition to the velocity distribution function shown in Fig. 5.4, we compute the velocity and angular correlation functions C v and C as a measure of the spatial range for which the velocity and orientation of a swimmer are coordinated with those of its neighbors [148]. The velocity and angular correlation functions are given by C v (r) = h( _ z i _ z j + _ z i _ z j )[rr ij ]i ij 2hj _ z i j 2 i i ; C (r) =hcos ( i j )[rr ij ]i ij ; (5.4) where r =jzj and r ij =jz i z j j. The angle bracketh:i represents the average over all possible pairs and the entire sampling time. These correlation functions are dened for r 2 since r < 2 lies inside the excluded area of the swimmer. These functions quantify how close the values of velocity and orientation of swimmers at dierent positions are, thus measuring the local order of the swimmers. To calcu- late the correlation functions numerically, we need to regularize the delta function properly. In this study, we discretize and regularize the sampling distance r of the correlation functions into multiple slabs, with each slab having a width of 0:2. That is to say, instead of measuring the correlation functions at single points, we compute the correlation functions based on the swimmers located inside the dierent slabs. 74 0 1000 2000 3000 4000 5000 0 0.2 0.4 0.6 0.8 1 t <P> r 0 0.1 0.2 0.3 0.4 0 0.2 0.4 0.6 0.8 Probability 2 3 4 5 6 7 8 9 10 0.6 0.65 0.7 0.75 0.8 0.85 0.9 r C v , C θ (a) (b) κ κ ν ν 1 2 <P> C v , C θ Probability 2 4 6 8 10 12 14 16 0 0.2 0.4 0.6 0.8 0 1000 2000 3000 0 0.2 0.4 0.6 0.8 1 t 0 0.1 0.2 0.3 0.4 0 0.1 0.2 0.3 0.4 (c) Φ A = 0.1 Φ A = 0.2 Φ A = 0.3 = 0.3, = 3 ν ν 1 2 Φ A = 0.1 Φ A = 0.2 Φ A = 0.3 = 0.3, = 0 Figure 5.6: Statistics associated with swirling (top row) and orientational order (bottom row): polar order parameterhPi, velocity and angular correlationsC v ,C , probability distribution of rotational activity parameter , all computed based on data corresponding to the shaded time range. Local correlations are observed for the swirling-type motion while global cor- relations are seen when orientational order is developed (Fig. 5.6(a)). The local correlations of the swirling motion are characterized by strong peaks located at po- sitions close to the swimmer surface, followed by a gradual decrease of correlation values with distance. This implies that the swimmers exhibiting swirling motions are aligned locally, not globally. In contrast to that, the swimmers exhibiting global orientational order have high correlation values over a long range of distance. In both cases, the values of correlation functions decrease when the area fraction A increases. 75 We also compute the polar order parameter hP (t)i = 1 N N X n=1 exp(i n (t)) (5.5) to assess the degree of global order in the swimmers population (Fig. 5.6(b)). The value ofhPi equals to 1 if all the swimmers are aligned in the same direction and is close to 1= p N if the population is randomly oriented. For the case of swimmers exhibiting swirling motions, the value ofhPi remains small over all time, whereas for the case of swimmers exhibiting global orientational order, the value ofhPi increases continuously and then levels o after long time. This indicates that the development of a global orientational order is process over a slow time scale, as opposed to the swirling and aggregation modes. This is because in the orientational order mode, the hydrodynamic interactions are dominated by the weaker gradient term in (5.1), as opposed to the velocity term, which explains why it develops at a much slower time scale. We also compute the long-time developedhPi, averaged over several samples, for dierent values of 2 and A while xing value of 1 , with the results summarized in Fig. 5.7(a). It is clear that the swimmers lose their global order whenj 2 j and A increases. The loss in global order at nite 2 lies on the fact that the orientation dynamics due to the gradient term in (5.1), which is responsible for the development of global order, is a higher order eect and is suppressed by the lower order orientation dynamics due to the alignment to local ow velocity. While the correlation functions and polar order parameter succeed in distin- guishing between the orientational order behavior and the swirling behavior, these 76 0.1 0.2 0.3 0 0.2 0.4 0.6 0.8 1 Φ A <P> 2 = 0.3 2 = 0 2 = 0.1 2 = -0.3 2 = -0.1 1 = 0.5 (a) 0 0.5 1 1.5 2 2.5 3 0 0.05 0.1 0.15 0.2 0.25 0.3 2 κ 1 = 0 Φ A = 0.3 Φ = 0.2 A Φ A = 0.1 (b) t → ∞ Figure 5.7: (a) The time-averaged mean and standard deviation of and (b) the sample-averaged mean and standard deviation of long-time developedhPi over var- ious initial conditions. statistical functions do not distinguish between the swirling behavior, for which the velocity and orientation are locally but not globally correlated, and the transitional behavior between the three global modes reported here, which is also characterized by locally-correlated velocity and orientation but no unambiguous swirling patterns. Therefore, we introduce a rotational activity parameter (t) = 1 N N X n=1 j _ n (t)j j _ z n (t)j (5.6) that measures the average change in orientation weighted by the traveling distance (Fig. 5.6(c)). The values of exhibit a continuous transition as 2 decreases as shown in Fig. 5.7(a). Fig. 5.7 also shows that the general dependence of and long-time developedhPi on 2 and A is insensitive to variations in initial con- ditions. We set threshold values for andhPi, and use the mean of the veloc- ity distribution function to distinguish between the three global modes: swirling, 77 −1 0 1 2 3 0.3 0.2 0.1 0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0.3 0.3 0.2 0.1 0 0.1 0.2 0.3 0.4 0.5 Φ A ν 1 ν 2 ν 2 Φ A ν 1 Swirling Orientational order Aggregation Figure 5.8: Phase diagrams showing the three dierent collective modes, swirling, orientational order and aggregation, as a function of 1 , 2 and A . Regions of transitional behavior for which no clear pattern are identied are left in white. orientational order, and aggregation. In particular, the swimmer populations are identied as exhibit swirling mode with time-averaged 0:15; when swimmer populations for a particular set of parameters have sampled averaged long-time de- velopedhPi t!1 0:5, the swimmers are identied as exhibiting the orientational order mode; and the aggregation mode is identied when the swimmer populations exhibit velocity distribution function with strong peak around zero velocity. These behaviors are mapped onto the 3D parameter space ( 1 , 2 , A ) in the phase dia- gram depicted in Fig. 5.8. Note that the orientational order mode arises in a narrow parameter range where 2 is small, which explains why the orientational order mode is usually not observed in experiments of collective behaviors of bacteria. To test the robustness of the global modes to physical randomness such as tum- bling eects of motile bacteria, we considered the eects of rotational Brownian noise on the global modes by adding a rotational diusion term Re D r ie in to 78 the orientation equation in (5.1). Here, D r is the rotational diusivity and fol- lows a normal probability distribution function with zero mean and unit variance. Note that, when modeling bacterial suspensions, it is customary not to account for translational diusion explicility [13, 69], since the translational diusion coe- cient for bacterial motion is much smaller than the rotational diusion coecient. Translational diusion is induced indirectly from rotational diusion through the coupling between orientation equation in (5.1) and the translation equation in (5.1). We set D r to have values of order 10 2 so that the corresponding P eclet number Pe = ( 1 + 2 )U A =(2D r ) is of order of 1. At this order of Pe numbers, the ad- dition of Brownian noise had negligible eects on the emergence of the swirling, orientational order and aggregation phase, causing only small perturbations to the \boundaries" of these global modes. 5.6 Conclusions To conclude, we note that, whereas our nite-sized system is fully-nonlinear and non-dilute, the ( 2 ; A )-slice of Fig. 5.8 may be interpreted as a rough analog to the 2D linear stability diagram obtained in the kinetic model of [13], see Fig. 1.5, where 1 was considered identically zero. The linear instabilities reported in the kinetic model do not reveal the nature of the emergent collective modes (swirling versus clustering) nor their physical implications. Here, we showed, for the rst time, that transitions from swirling to clustering and aggregation, and vice-versa, can be induced by an interplay between the agellar activity and the hydrodynamic 79 interactions. These results, albeit in the context of a simplied model, may have profound implications on understanding intermediate stages where agellar activity is decreased in response to environmental conditions such as, overcrowding or nutri- ent depletion. Our ndings suggest that the interplay between hydrodynamics and agellar activity may serve as a physical mechanism for aggregating the bacteria, thus altering their convective ow which may, in turn, alter microbial processes such as transport of oxygen and nutrients. 80 Chapter 6 Circularly conned microswimmers exhibit multiple global patterns 6.1 Introduction The eects of geometric connement on the dynamics of individual swimmers and swimmer populations are becoming the center of attention of several experimental and theoretical studies. At the individual level, the eects of wall boundaries on cell locomotion have been intensively studied. For example, bacteria with helical ag- ella such as E. Coli, were shown to change their swimming trajectory from straight to circular near a surface [66]. Besides, wall-induced hydrodynamic interaction can lead to attraction or repulsion of a swimmer from a solid planar surface, depending on the whether the swimmer is a pusher or puller. This hydrodynamic induced attraction provide an explanation for accumulation of swimming cells near surfaces, which was observed in many biological experiments [9]. However, at the population level, little is known about the eects of geometric connement on the emergent 81 global patterns in the self-propelled systems. It is only recently that systematic experiments were performed to investigate the eects of circular connement on the emergent behavior in the self-propelled system. For examples, recent experi- ments on conned bacterial suspensions within a cylindrical droplet show that, if the radius of connement is below a critical value, the bacteria can spontaneously or- ganize themselves into a spiral vortex encircled by a counter-rotating cell boundary layer [142], see Fig. 1.2(a). This remarkable self-organization is attributed to the in- terplay between hydrodynamic and boundary interactions [77]. Other experiments on conned bacterial suspensions in a spherical water droplet report similar results, namely, a global vortical motion at relatively high bacterial concentration, albeit with no xed axis of rotation. They also report that, at relatively small concen- trations, the bacteria tend to aggregate at the inner boundary of the droplet [136]. This emergence of coherent vortical motion is not limited to bacterial suspensions. Recent experiments on electrically powered motile colloids in a Hele-Shaw cell with circular connement also reported the formation of vortex state by the colloids, albeit the colloidal vortex is characterized by a low density center instead of a counter-rotating swimmer layer [11], see Fig. 1.2(b). The experimental evidence in [11, 136, 142] suggests that geometric connement aects the global order of microswimmers. While motivated by these studies, our goal here is not to theoretically reproduce their results. Rather, it is to understand how the interplay of both Hele-Shaw and circular connement in uence the emer- gent collective dynamics in the context of an idealized model. We start by deriving 82 a closed-form expression for the velocity eld created by N such swimmers in a circularly-conned domain. We then investigate numerically the emergent behav- ior aorded by a population of microswimmers subject to Hele-Shaw and circular connement as a function of the connement radius and agellar activity. We show that decreasing agellar activity induces a hydrodynamically-triggered transition in conned microswimmers from swirling to global circulation (vortex) to bound- ary aggregation and clustering. These results highlight that the complex interplay between connement, agellar activity and hydrodynamic ows in concentrated suspensions of microswimmers could lead to a plethora of global patterns that are dicult to predict from geometric consideration alone. 6.2 Hydrodynamic model for conned microswim- mers in circular domain Following the result in Chapter 5, the dynamics of a population ofN microswimmers in Hele-Shaw connement can be expressed in concise complex notation as _ z n =Ue in +w(z n ) +V n ; _ n = Re wie in ; (6.1) Here,z n denotes the swimmer's position in the complex planez =x+iy (i = p 1) and n its orientation, n = 1;:::;N. The operator Re[] takes the real part of the expression in bracket whereas the overline () denotes the complex conjugate 83 α x y h R R R / z z α π−α n n 2 x y (a) (b) Figure 6.1: Microswimmers in Hele-Shaw and circular connement: (a) schematic of the set-up and the dipolar far-eld ows created by a swimmer in Hele-Shaw connement (inset); (b) Eect of the circular connement is accounted for using an image system outside the circular domain as dictated by the Milne-Thompson circle theorem. of the expression in parenthesis. Each swimmer has a self-propelled speed U and is advected by the local ow velocity w. It is understood that all variables and parameters are non-dimensionalized based on characteristic length and time scales dictated by the individual swimmers. The value of re ects the degree of agellar activity (or geometric polarity) and could decrease from a positive value ( > 0) for vigorously-beating agella (or large-tail particle) to a negative value ( < 0) for weakly-beating agella (or large-head particle), thus causing the swimmer to switch from reorienting with the local ow velocity to reorienting in the direction opposite to the local ow. For simplicity, the reorientation of swimmers in response to the local ow gradient is neglected, since its eect is small when the absolute value of is suciently greater than zero as shown in Fig. 5.8 in Chapter 5. To take into account the eect of circular connement on the hydrodynamic signature of the microswimmers, we apply the Milne-Thomson circle theorem, see, 84 for example, [4]. This theorem states that for any irrotational ow with an associ- ated complex potential functionF (z), analytic in the considered domain, a circular boundary of radius R can be constructed at the origin using the modied complex potential function G(z) =F (z) +F (R 2 =z): (6.2) The second term in the expression enforces the imaginary part of the modied complex potential functionG(z) to be zero atjzj =R, thus creating an impenetrable circular boundary of radius R. In Hele-Shaw connement, the leading-order ow disturbance created by each swimmer is that of a potential source dipole, see Fig. 6.1(a) and [13, 69, 130]. The complex potential function for a system of N swimmers is therefore given by F (z) = N X j=1 e i j zz j ; (6.3) where the dipole strength =a 2 U is proportional to the swimming speedU and the square of the eective radiusa of the swimmer, which are here normalized to 1. The complex conjugate velocity at each swimmer's position, excluding the self-induced ow disturbance, is given by w(z n ) = dG dz z=zn e in (zz n ) 2 : (6.4) 85 By virtue of (6.2) and (6.3), Eq. (6.4) becomes w(z n ) = N X j6=n e i j (z n z j ) 2 + N X j=1 (R 2 =z 2 j )e i( j ) (z n R 2 =z j ) 2 : (6.5) That is, constructing the ow eld created by a swimmer located at z n and ori- entation of n in a circular boundary of radius R amounts to placing an image swimmer at locationR 2 =z n , with orientation n and dipole strength(R 2 =z 2 n ). A schematic of the image system is illustrated in Fig. 7.1(b). We introduce a collision avoidance mechanism among swimmers and between a swimmer and the circular boundary based on the rapidly-decaying repulsive part of the Leonard- Jones potential V n = N X j6=n f 0 " 2a jz n z j j 13 z n z j jz n z j j + a jdj 13 z n d jz n dj # : (6.6) Here, f 0 is a scaling parameter that characterizes the repulsion strength, d is the shortest vector joining the swimmers and the circular wall (parallel to the perpen- dicular bisector of the circle). The rst term ensures collision avoidance among swimmers whereas the second term guarantees that no swimmer penetrates the circular boundary. These near-eld interactions decay rapidly outside a small ex- cluded area centered around z n . Their rapid decay ensures that the order of the far-eld hydrodynamic interactions is preserved. Equations (6.1), (6.5) and (6.6) form a closed system of equations forN conned swimmers in a circular domain of 86 radius R that take into account the swimmer-generated ows and level of agellar activity. Note that the model presented here diers signicantly from the model in [137] which account for the interactions between the microswimmers and environment using an eective potential force without consideration of hydrodynamic eects. While closer to the model in [77], the two models are distinct in the underlying physics and technical implementation. The leading order ow disturbance created by the swimmers are not of the same order: here, swimmers create a source dipolar eld that decays as 1=r 2 with distance r from the swimmer. In [77], solving Stokes ow in 2D domains, a swimmer's force dipolar ow decays as 1=r. On the technical side, here we derived closed-form expressions that account for the eect of the circular connement on the hydrodynamic interaction between the micro-swimmers, as opposed to the approximate numerical approach in [77]. The swimmers here are subject to friction from the upper/lower walls while these eects are neglected in [77]. Steric interactions are added here to avoid direct collision but the orientation response in equation (6.1) is purely hydrodynamic. In [77], steric interactions in the form of Gay-Berne potential were used to account for both collision avoidance and orientation of the swimmers at the circular boundary. 87 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 (a) (b) (c) −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 −50 −25 0 25 50 Figure 6.2: Snapshots for the emergence of three global modes: (a) chaotic swirling ( = 3); (b) stable circulation ( = 0:5); (c) boundary aggregation ( =0:5). The corresponding ow elds are depicted at the bottom row. For all three simulations, N = 750 and R = 50. 6.3 Emergence behaviors: chaotic swirling, global vortex circulation, boundary aggregation We examine the emergent global modes of populations ofN swimmers as a function of the connement radius R and agellar activity . We normalize U and to 1 and set = 0:5. We x the area fraction of the swimmers = Na 2 =R 2 to 0.3, that is, we vary the number of swimmers N in the domain as R changes. We focus on initial conditions for which the swimmers are initially randomly oriented 88 but spatially homogeneous. We perform Monte-Carlo type simulations in the sense that for each set of parameters, we run multiple trials corresponding to dierent sets of initial conditions, each set taken from a uniform probability distribution function. We observe that, as we vary the level of agellar activity for xed radiusR, the swimmers organize into three distinct types of global modes: chaotic swirling, stable circulation and boundary aggregation. Snapshots for representative simulations are shown in Fig. 6.2, together with their corresponding instantaneous velocity elds. The emergence of chaotic swirling mode was observed in our previous study on conned swimmers in a doubly-periodic domain in Chapter 4, and was explained on the ground that swimmers have vigorously-beating agella ( positive) tend to tail-gate other swimmers ahead and form chain-like structures. This results in vortices and swirls in the velocity eld, see Fig. 6.2(a) left, that is reminiscent of the bacterial turbulence observed in numerous experiments [21, 34, 36, 81, 114, 115]. A surprising emergent mode arises when we reduce to smaller positive values. In this case, the swimmers spontaneously organize into stable circulation, charac- terized by a single vortex with swimmers concentrating at the circular boundary, see Fig. 6.2(b). The circulation occurs with equal probability in clockwise and anti- clockwise direction. Similar phenomenon was observed in the recent experiments and simulations of bacterial suspensions in conned circular domains, [77, 142] and self-propelled colloids in a circular Hele-Shaw cell [11]. However, in contrast to the single vortex state presented here, [77, 142] report a pair of counter-rotating 89 vortices, constituting a single vortex of swimmers encircled by a counter-rotating layer of swimmers. The model in [77] was designed to explain the physical origin of the counter-rotating layer of swimmers near the circular boundary. The simpler model presented here is able to capture the leading-order global behavior observed in [77, 142] which is that of stable circulation. It is worth noting that our result of single vortex of microswimmers with a low density center agrees with the vortex structure observed in the experiments of motile colloid rollers in a Hele-Shaw type circular geometry [11], which perhaps is a closer set up to our study. In [11], the colloids are subject to both hydrodynamic and electrostatic interaction. Albeit the interaction among the motile colloids is not purely hydrodynamic, hydrodynamic interaction should play an important role in the formation of the vortex structure. As decreases to a negative value, indicating swimmers with even weaker ag- ellar beating, the swimmers aggregate at the circular boundary or form clusters inside the domain, see Fig. 6.2(c). The formation of `inner' clusters was observed in our previous study on conned swimmers in a doubly-periodic domain in Chapter 5. It was explained based on the fact that swimmers with weakly-beating agella tend to reorient in the opposite direction to the local ow created by nearby swim- mers and, thus, tend to aggregate and cluster. These clusters have the shape of an asterisk and create a source-like convective ow that attracts other swimmers to the cluster, potentially causing it to develop into large cluster. The aggregation inside the circular boundary can be attributed to a similar phenomenon. However, 90 the aggregation at the boundary is due to both steric interactions as well as to the swimmers tending to cluster with their image system outside the circular domain. 6.4 Statistics of global modes We use a number of statistical measures to assess the observed global structures. We focus on the statistics of the circulation mode. The statistics of the other two modes (swirling and aggregation) have been described in details in our previous study in Chapter 5 and are omitted here for brevity. Following [142, 77], we introduce the vortex-order parameter = P j v t;j = P j j _ z j j 2= 1 2= ; (6.7) where v t;j = _ z j t j + _ z j t j =2 is the tangential speed of each swimmer and t j is the vector parallel to the tangent of the circle at the location of the swimmer. one gets = 1 for steady circulation and = 0 for disordered state. The change in with time is presented in Fig. 6.3(a) for the case shown in Fig. 6.2(b). The value of levels o when the swimmers develop into a stable circulation. Fig. 6.3 also depicts the mean tangential speed and mean polarization eld for the same case. We compute the mean tangential speed and the probability distribution of the swimmers as a function of the distance r = p x 2 +y 2 from the origin, see Fig. 6.3(b). The tangential speed increases withr, indicating that the swimmers have a higher vortex order closer to the circular boundary. The probability distribution also increases 91 with r, re ecting that the vortex forms near the inner boundary. The polarization eld shows the average orientation of the swimmers at dierent location in the vortex state and is computed by interpolating the orientation of the swimmers to a polar grid. The center of the polarization eld, where no clear orientation can be determined, is left out in gray color. The formation of vortex structure can be seen clearly from the polarization eld as shown in Fig. 6.3(c). We next examine the eect of the interplay between agellar activity and connement radiusR on the developed circulation as measured by the vortex-order parameter . Fig. 6.4(a) shows the mean long-time for several initial conditions and for dierent values of andR. AsR increases, the value of decreases almost linearly when R < 60 and the rate of decrease in increases with . Large error bars are observed in the cases of R = 60 and R = 70 for = 1:5, which mark the presence of an instability that prevents the swimmers from forming stable cir- culation. This indicates that the swimmers can only form stable circulation when the radius of the circular boundary is below a critical value whose precise value de- pends on the agellar activity. The existence of critical radius for stable circulation is also observed in the aforementioned experiments of bacterial suspensions within a circular droplet [77, 142]. In their case, the value of rst increases withR, then levels to a constant value for a range ofR beforeR exceeds a critical value for which suddenly drops. In our case, is monotonically decreasing as a function of R before R reaches a critical value for which no stable circulation is observed. 92 -50 0 -25 50 25 -50 -25 0 25 π/2 π 3π/2 2π 5π/2 50 0 200 400 600 800 1000 0 0.2 0.4 0.6 0.8 1 t Φ (c) (b) (a) 0 x y 0.6 0.8 1 1.2 1.4 0 0.05 0.1 0.15 0.2 v t Probability 10 20 30 40 50 r Figure 6.3: Statistical data for the stable circulation mode shown in Fig. 6.2(b): (a) vortex order parameter ; (b) mean tangential speed and probability distribution of the swimmers; (c) mean polarization eld. The mean polarization eld and tangential speed are computed based on the averaging of data over the shaded time range, which a stable vortex order has been developed. 93 10 20 30 40 50 60 70 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Φ t→∞ (a) R 10 20 30 40 50 60 70 −1 −0.5 0 0.5 1 1.5 2 2.5 3 R ν (b) chaotic swirling stable circulation boundary aggregation ν = 0.5 ν = 1 ν = 1.5 Figure 6.4: (a) The time-averaged mean and standard deviation of long-time devel- oped over various initial conditions. The standard deviations are denoted by the error bars. (b) Phase diagram showing the three dierent global modes: chaotic swirling, stable circulation and boundary aggregation, as a function of and R. We set threshold value of t!1 = 0:5 to distinguish between the stable circula- tion and chaotic swirling modes. We then map three global modes: chaotic swirling, stable circulation and boundary aggregation onto the 2D parameter space (, R) depicted in Fig. 6.4(b). This plot summarizes the main ndings of this work and identies the critical values of circular connementR that mark the transition from chaotic swirling to stable circulation for dierent values of agellar activity . 94 6.5 Conclusions Our minimal model for microswimmers in Hele-Shaw and circular connement ac- counts for swimmer- uid interactions, swimmer-generated ows, and level of agel- lar activity and captures qualitatively the dynamics seen in experiments on conned bacterial suspensions. In particular, we showed the following: (i) For vigorously-beating agella, the global swirling pattern is reminiscent to turbulence-like patterns observed in many experiments on bacterial suspensions. (ii) As the agellar activity or the connement radius decreases, the microswim- mers stabilize into a global vortex with swimmers' concentration increasing towards the circular boundary. (iii) As the agellar activity decreases further so that the swimmers begin to reorient in the opposite direction to the local ow, the swimmers spontaneously cluster and aggregate inside and at the circular boundary. This interesting global phenomenon is reminiscent to the recently reported experimental results on bacte- rial cells aggregating at the boundary of a spherical droplet [136]. These transitions in the global order of conned microswimmers lend a note of caution when drawing conclusions on emergent patterns, especially in experiments on living organisms and bacterial suspensions. Experiments are often conducted under specic environmental conditions and parameter values. Thus, the observed global structures may only be valid under those conditions. For example, while one may conclude from [142] that suitably designed geometric boundaries provide 95 a means for stabilizing and controlling order in active microbial systems, this con- clusion overlooks the fact that time-dependent bacterial properties, such as level of agellar activity, may in uence and alter the global order. Future extensions of this work will account for the eects of the shape asymme- try on the steric interactions as in [140] in order to discern the interplay between hydrodynamic and steric asymmetry on the emergent dynamics. 96 Chapter 7 Density shock waves in conned microswimmers 7.1 Introduction The transport of micro-particles in ow channels occurs in a variety of industrial and natural processes. Examples include ltration [48], colloid deposition [30], droplet-based micro uidics [8], blood cells in micro ows [94], sperm cells in the fallopian tube [99], and pathogens in the urinary tract [86]. The ability to control the density distribution and group velocity of micro-particles in ow would therefore have numerous applications in physics and biology. This chapter analyzes, via discrete particle simulations and macroscopic models, the emergent patterns in populations of motile particles driven by an external ow in a rectangular micro- channel. The dynamics of single and many-swimmer systems are typically studied in un- conned settings [36, 79, 105]. The eects of connement to narrow ow channels 97 are less well understood; see, e.g., [9, 40, 42, 46, 150] and references therein. Recent experiments on driven droplets [6, 7, 8, 32] and self-propelled colloids [11, 12] in quasi two-dimensional (Hele-Shaw type) channels show the emergence of traveling density waves, including density shocks at the wave front, see Fig. 1.4. Similar observations are reported in simulations of conned self-propelled particles [69], al- beit with density shocks forming at the back of the wave.The primary mechanism responsible for the emergence of these density waves is attributed to hydrodynamic interactions (HI) among the conned particles [12, 69]. Connement in Hele-Shaw type geometries induces a distinct hydrodynamic signature; the particle far-eld dis- turbance is that of a source dipole irrespective of the details of the particle transport mechanism, driven or self-propelled, pusher or puller [6, 8, 7, 13, 33]. Additional connement to a narrow channel does not alter the nature of hydrodynamic in- teractions but imposes impenetrability conditions on the lateral boundaries of the channel. These eects are properly accounted for in the model used in [69], and, thus, do not explain the discrepancy in the location of the shock formation between the simulations [69] and experiments [6, 12]. Here, we show that the background ow and geometric connement coupled with the HI among the swimmers and their motility properties lead to a rich variety of emergent behaviors, consistent with ex- perimental evidence but unobserved in simulations before. Most notably, we show the emergence of density shocks that transition from `subsonic' to `supersonic' as the intensity of the background ow increases. We conclude by demonstrating that 98 α x y n y n W −αn −α −α W− − W W y − − n n y y 2W h W V L (a) (c) (b) U 2W V α + Figure 7.1: (a) Schematic of microswimmers in a Hele-Shaw rectangular channel. (b) Dipolar far-eld ows induced by the self-propelled motion and the background ow. (c) Sidewall connement is accounted for using an innite image system. this richness can be exploited to control the density distribution and collective speed of the microswimmers. 7.2 Hydrodynamic model for conned microswim- mers in rectangular channel For concreteness, we brie y review the model of N microswimmers in Hele-Shaw connement. Namely, one has _ z n =Ue in +w(z n ) +v n ; _ n = Re wie in : (7.1) 99 Here,z n denotes the swimmer's position in the complex planez =x+iy (i = p 1) and n its orientation,n = 1;:::;N. The overline () and the operator Re[] denote the complex conjugate and the real part, respectively. Each swimmer has a self- propelled speed U and is advected by the local ow velocity w(z). This model is applicable when the swimmer size, which we denote by 2a, is of the same order as the Hele-Shaw heighth (Fig. 7.1), that is, whenh=2aO(1). In (3.17), all variables are non-dimensional using the characteristic length and time scales dictated by the size a and speed U of the individual swimmer. and are translational and rotational mobility coecients whose values depend on the geometric properties of the swimmer and its agellar activity: is in the range 0 < < 1 while is either positive or negative depending on the head-tail asymmetry [13] or agellar activity [130]. We use a collision avoidance mechanism v n based on the repulsive part of the Leonard-Jones potential v n =f 0 " N X j6=n 2R jz n z j j 13 z n z j jz n z j j i R W=2jy n j 13 y n jy n j # : (7.2) f 0 is the scalar parameter that characterizes the strength of repulsion, and we set f 0 = 1. These near-eld interactions decay rapidly outside a small excluded area centered around z n , thus ensuring that the order of the far-eld HIs is preserved. To close the model in (7.1), we need to evaluate the velocity eldw(z) created by N conned swimmers in uniform background ow, (Fig. 7.1(a)). The self-propelled motion of a swimmer z j induces a dipolar velocity eld w s (zz j ) =a 2 Ue i j =(z z j ) 2 which points in the swimming direction, whereas the background ow induces 100 another dipolar eldw b (zz j ) =a 2 (1)V=(zz j ) 2 which points in the direction opposite to the background ow, irrespective to the orientation of the swimmer (see Fig. 7.1(b) and [33]). Here,a andU are the eective radius and self-propelled speed of the swimmer, normalized to 1, andV is the velocity of the background ow. This model is similar to the one used in [69]. However, the fact that w b is proportional to (1) was overlooked in [69], thus failing to correctly capture the dipolar eld induced by the background ow. The velocity eld w(z) created by N swimmers in uniform background ow can be written as a superposition of w s and w b . The impenetrability conditions at z =iW=2 due to the sidewall connement are accounted for using the method of images. Basically, to each swimmer at z n = x n + iy n correspond two periodic lattices, of period 2iW , of image dipoles located at z n =x n + iy n and z n =x n + i(Wy n ), respectively (Fig. 7.1(c)). To calculate the velocity eldw due toN swimmers and their image system, one has to evaluate innite sums of terms that decay as 1=z 2 . An approximate numerical solution is given in [69]. Here, we derive an exact analytical solution in terms of hyperbolic functions. The complex potential of a source dipole induced by a swimmer at the origin in an unbounded domain is given by F s =a 2 Ue i =z. It follows that in a singly periodic domain with period 2iW , the corresponding complex potential is given by, F periodic =a 2 Ue i 1 z + 1 X n=1 1 z 2inW + 1 z + 2inW ! : (7.3) 101 Recall the series representation cotz = 1 z + P 1 n=1 1 zn + 1 z+n . Equation (7.3) can then be rewritten in a closed-form expression, F periodic = a 2 U 2W e i coth z 2W : (7.4) Now we consider N swimmers with the image system depicted in Fig. 7.1(c). The location of the n-th swimmer is denoted by z n =x n + iy n . Using (7.4), and taking into account the eect of uniform background ow V and translational mobility coecient , we obtain the complex potential for the image system, F =Vz a 2 2W ( N X j=1 [Ue i j (1)V ] coth h 2W (zz j ) i + [Ue i j (1)V ] tanh h 2W (zz j ) i ) : (7.5) where j =Ue i j (1)V anda 2 j indicates the strength and orientation of the dipole induced at z j . The complex conjugate velocity eld due to the swimmers can then be evaluated through w = dF=dz. Excluding the self-induced ow disturbance, we obtain the complex conjugate velocity eld at the location of the n-th swimmer, w(z n ) =V + a 2W 2 ( N X j6=n j csch 2 h 2W (z n z j ) i N X j=1 j sech 2 h 2W (z n z j ) i ) (7.6) where j = Ue i j (1)V and a 2 j indicates the strength and orientation of the dipole induced at z j . When the swimmers are subject to a suciently strong 102 background ow, they quickly reorient and align in either the same (for > 0) or opposite (for < 0) direction to the background ow. Basically, the background ow suppresses the orientation dynamics of the swimmers as noted in [69]. All j reach, in nite time, a developed value d = pU (1)V , where p = sgn() and sgn is the signup function. It turns out that d , which we refer to as the eective polarization parameter, is a key parameter for understanding the emergent structure of the swimmers. 7.3 Speed of sound To understand the propagation of density shock waves of microswimmers in the conned system, we need to rst identify the speed of sound in the conned medium. The speed of sound in a medium is dened as the speed of propagation of a small amplitude density perturbation in that medium. What is the speed of sound in our medium? To answer this question, we rst have to identify the medium. Here, we distinguish three cases: (a) a ow channel where the population of microswimmers is very dilute such that the hydrodynamic interactions (HIs) among the swimmers are negligible, (b) a ow channel where a dense plug of microswimmers lls part of the channel, which is the main topic of our study, and (c) a densely lled channel, see Fig. 7.2. Cases (b) and (c) are considered in [6, 8] and in [32], respectively, albeit for driven droplets. In [6, 8], compression shock waves develop at the front of the droplet population. These shock waves are called supersonic because the speed of the density wave 103 formed by the droplets is faster than the speed of sound in the medium ahead. Following the same terminology, we refer to the compression shock waves observed in our model as subsonic when the wave speed is slower than the speed of sound in the medium ahead and as supersonic when it is faster than the speed of sound in the medium ahead. V (a) dilute channel (b) dense plug in dilute channel (c) dense channel no shock no shock shock formation c = U + μV c = U + μV c = U + μV + .... Figure 7.2: The gradient in the density distribution of the initial swimmer popula- tion is essential for the formation of the compression shock waves. Here, we outline an analytical approach for obtaining the speed of sound from a continuum density model, which is derived by following a standard homogenization procedure. We start by considering the 2D continuum density r (r) with r =x^ x+y^ y and ^ x; ^ y are unit vectors in Cartesian coordinates. The continuum equation that governs the evolution of the density function r (r) is given by @ r @t +r r nh (Up +V )^ x + Z w c (r; r 0 ) r (r 0 )dr 0 i r (r) o = 0: (7.7) Using polar coordinates representation, one has x = r cos, y = r sin. For a channel of innite width, the dipolar velocity eld is given by w c (r; r 0 ) = d (rr 0 ) 2 (cos 2^ x + sin 2^ y): (7.8) 104 Any dipolar velocity contributions within the excluded areajr r 0 j < 2a are re- moved from the integral in (7.7). We start from a constant density o and perturb and linearize (7.7) around o , namely, we set r = o + ~ , (7.7) then becomes @ ~ @t + @ @x [(Up +V )~ ] +r r (I 0 +I 1 ) = 0; (7.9) where I 0 = o ~ Z w c (r; r 0 )dr 0 ; I 1 = o Z w c (r; r 0 )~ (r 0 )dr 0 : (7.10) In the dilute case where HIs are negligible,I 0 andI 1 are identically zero and the speed of sound of the medium isc =Up +V . This is the case shown in Fig. 7.2(a) and also the medium ahead of the dense plug of swimmers in Fig. 7.2(b). For completeness, we also report the speed of sound for the case shown in Fig. 7.2(c) by analytically computing I 0 and I 1 . I 0 = o ~ d Z 2 2 Z L 1 cos 2a + Z 3 2 2 Z L 2 cos 2a cos 2 r 2 rdrd = o ~ d : (7.11) Note that the solution is independent of L 1 and L 2 , indicating that it is valid for any arbitrary slab of swimmers in the channel. We now consider the integral I 1 . We take the divergence of I 1 , r r I 1 = o r r Z w c (r; r 0 )~ (r 0 )dr 0 : (7.12) 105 Applying divergence theorem to (7.12), we get r r I 1 = o 2a Z 2 0 w c (r; r 2a^ r 0 ) ^ r 0 ~ (r 2a^ r 0 )d; (7.13) where ^ r 0 = cos^ x + sin^ y. Equation (7.13) can be simplied as r r I 1 = o d 2a Z 2 0 ~ (r 2a^ r 0 ) cosd: (7.14) Consider the perturbation in the form of plane wave, ~ (r) =Ae j(kr!t) =Ae j(kx+ly!t) . We focus on 1D perturbation and set l = 0. Direct substitution into (7.14) yields r r I 1 = A o d 2a e j(kx!t) Z 2 0 e j2ak cos cosd: (7.15) The resulting integral in (7.15) can be evaluated using Bessel function of the rst kind [1], that is, r r I 1 = A o d j a J 1 (2ak)e j(kx!t) : (7.16) We obtain the dispersion relation of the sound by substituting (7.11) and (7.16) into (7.9), ! =k(Up +V o d ) o d a J 1 (2ak): (7.17) 106 The speed of sound in a homogenous medium is therefore given by For non-zero wavenumber k, the speed of sound is given by c = ! k =Up +V o d o d ak J 1 (2ak); k> 0: (7.18) In the long wave limit, lim k!0 c =Up +V o d [1 +J 0 (0)J 2 (0)]: (7.19) Here,J 0 ,J 1 andJ 2 are Bessel functions of the rst kind. Note that a localized per- turbation, as opposed to the periodic perturbation, travels with the group velocity c g = @! @k =Up +V o d [1 +J 0 (2ak)J 2 (2ak)]: (7.20) The phase velocity and group velocity of the perturbation match to the same value at the long wave limit. By way of validation, we also calculate numerically the speed of sound in a channel with nite width using a linearized version of the 1D continuum equation model described in the main text of the letter. It shows good agreement with the analytical result we presented above. Fig. 7.3 depicts the speed of soundc as a function of the wavenumberk in the moving frame of reference with velocity Up +V . 107 0 1 2 3 4 5 6 0.8 1 1.2 1.4 1.6 k W = ∞ (analytical) W = 60 (computational) c / ρ o 0.6 Figure 7.3: Speed of sound in the moving frame of reference with velocityUp +V . Parameters used: d = 1, = 0:5. 7.4 Subsonic to supersonic transition in density shock waves We examine the emergent behavior of populations of N swimmers that initially ll the whole channel width and a segment ` of the channel length and are homoge- neously placed and randomly orientated. Periodic boundary conditions are imposed at the two ends of the channel. We set = 0:5. The swimmers self-organize into a shock wave structure where, for > 0, the location of the compression wave depends on the strength of the background ow V , which in turn dictates the `speed of sound' of the medium given byc =U +V . For a weak background owV <c, the eective polarization d of the swimmers is dominated by their self-propelled motion so that d > 0. The emergent structure is characterized by a compression wave at the back and a rarefaction or expansion 108 wave at the front (Fig. 7.4(a)). The average population speed is smaller than c, see Fig. 7.5, thus this shock wave can be termed as `subsonic'. If we tune the background ow properly such that the polarizations induced by the self-propelled motion and background ow cancel out such that d = 0, the shock wave is suppressed and the swimmers evolve as a uniform cluster traveling at the speed V (Fig. 7.4(b)). When the background ow is strong V > c, it dominates the swimmers' eective polarization d < 0; the density shock forms at the front rather than at the back (Fig. 7.4(c)), and the average population speed is greater than c, see Fig. 7.5, thus the term `supersonic'. Note that for < 0, the swimmers align in the opposite direction to the background ow. The developed polarization parameter d is al- ways negative, and thus the swimmers develop a shock wave structure similar to Fig. 7.4(c) that could travels upstream or downstream depending on the value of V . The density shock shown in Fig. 7.4(c) is reminiscent to the experimental results of [8] for driven droplets in a Hele-Shaw channel and of the propagating density bands shown in [12] for self-propelled colloids. The propagating density fronts in the latter system and their similarity to our results are surprising; however, the agreement with [8] is consistent with the intuition that active swimmers behave like passive particles when subject to strong background ow. 109 −30 0 30 −30 0 30 0 200 400 600 −30 0 30 x y 0 0.4 0 0.4 0 0.4 0 200 400 600 x ρ shock at back shock at front shock suppressed V = 1 V = 2 V = 3 σ d = 0.5 σ d = 0 σ d = -0.5 ρ ρ (a) (b) y y (c) = 1.5 = 2 = 2.5 c c c Figure 7.4: Density shock waves in conned microswimmers. Left: (a) for V < c, compression shock wave forms at the back; (b) for V =c, shock is suppressed; (c) forV >c, compression wave forms at the front. The `speed of sound' of the medium ahead is c = U +V . Right: Comparison of simulation results between discrete particle and continuum model. The density is the local area fraction. Simulation results of the discrete particle model are colored in blue. Simulation results of the quasi 1D continuum model are superimposed in red. Parameter values are: W = 60, A = 0:3, = 1, snapshots taken at t = 1000. 0 0.2 0.4 0.6 0.8 1 Probability 1 1.5 2 2.5 3 swimming speed (U + μ w) shock at back shock at front shock suppressed Figure 7.5: Comparison of the speed of the density shock wave and the speed of sound for the case shown in Fig. 7.4. The probability distribution of swimming speeds of the microswimmers in the discrete simulations are shown in solid lines, and the speed of sound of the medium ahead is indicated in dashed lines. 110 shock at front shock at back σ > 0: d σ < 0: d V < c V > c Figure 7.6: 1D model does not capture correctly the shock wave formation. It results in compression shock at the front for V < c and at the back for V > c, as opposed to the 2D system which exhibits compression at the back for V < c and front for V >c (Fig. 7.4). The density shock can also be visualized by the local area fraction (x) of the microswimmers, depicted in the right panel of Fig. 7.4. The density (x) from the discrete particle simulations is computed as in [8], namely, we calculate the area occupied by the microswimmers in discretized slabs, with each slab having width of 2a. For each slab, the local area fraction is given by (x) = N x a 2 =(2aW ), wherex is sampled from discrete points at the centerline of the slabs and N x is the corresponding number of swimmers inside the considered slab. We now consider the simple picture of three swimmers perfectly aligned parallel to thex-axis. One can easily see that forV <c and d > 0, HI push the swimmers forward such that the middle swimmer is always faster leading to compression at the front, whereas for V > c and d < 0, HI hinder the swimming motion and the middle swimmer is always slower leading to compression at the back (Fig. 7.6). This argument can be generalized toN swimmers. Thus, in 1D systems, the density shock forms at the front of the swimmers when the background ow is weak and at the back when the background ow is strong. This prediction, illustrated in Fig. 7.6, agrees with the density shocks seen in the 1D experiments of driven micro uidic 111 droplets [17], and is exactly opposite to the 2D experiments [8] and simulations presented in Fig. 7.4. The contradiction reveals the importance of 2D HI and the complexity of the mechanism responsible for the formation of density shock waves in populations of microswimmers in conned 2D channels. Density shocks emerge as a result of the interplay between 2D hydrodynamics, sidewall connement, and the initial distribution of swimmers into a segment ` of the full channel. ForW =1 and swimmers homogeneously lling the whole space, no shock develops due to the radial symmetry of the dipolar eld at each swimmer. When sidewall connement is introduced, it serves as a mechanism to break the radial symmetry of the dipolar eld and creates a non-zero average hydrodynamic interaction acting on the swimmers. The initial uniform distribution of swimmers into a segment` of the full channel length creates a spatial gradient in the HI in the x-direction, which contributes to triggering the instability of the initially uniform swimmer distribution. We numerically investigate the eect of 2D interactions on the shock wave for- mation by xing = 1 and V = 1 (i.e., d = 0:5) and varying the channel width W , in unit ofa, and the initial local area fraction A =a 2 N=(W`). To obtain the desired value of A , we set ` = 100 and vary N. For each set of parameters, we run multiple trials corresponding to dierent sets of initial conditions taken from 112 Φ A W x 0 100 200 300 ρ(x) 0 0.1 0.2 0.3 t 0 500 1000 1500 2000 S 0 0.2 0.4 0.6 0.8 (a) (b) (c) 0.1 0.15 0.2 0.25 0.3 20 40 60 80 100 −0.2 0 0.2 0.4 0.6 < S > t→∞ Figure 7.7: (a) Shock order parameter is dened by comparing the actual density prole to a triangular density prole, shown here at t = 2000 for W = 80 and A = 0:3. (b) Evolution of the shock order parameter for the case in (a). (c) Mean value of the shock order parameter as a function of channel widthW and swimmers' area fraction A , computed base on the shaded time range and averaged over all sample simulations. The dashed line shows the contour line ofhSi t!1 = 0:4, above which shock waves are unambiguously observed. a uniform probability distribution function (Monte-Carlo type simulations). We capture the emergence of the density shock wave using the shock order parameter S(t) = " 1 2 X x j(x;t) t (x;t)j # = X x t (x;t); (7.21) where(x) is the local density of the swimmers and t (x) is a triangular t of the density prole(x) (Fig. 7.7(a)). is a signed parameter with value determined by the slope of the inclined line of the tting triangle: = 1 corresponds to a negative 113 slope and =1 corresponds to a positive slope. S = 1 corresponds to a perfectly triangular subsonic shock, whereas S =1 corresponds to a perfectly triangular supersonic shock, and S = 0 to a rectangular homogenous prole. An example of the evolution of S is depicted in Fig. 7.7(b). The dependence of the mean value hSi t!1 averaged over all simulations on the parameter space (W; A ) is depicted in Fig. 7.7(c), which indicates the values of W and A necessary for the shock waves to develop. In particular, it shows a continuous increase inhSi t!1 as W and A increase which then levels o for largeW and A . This plot provides two important pieces of information. First, it serves as a strong supporting evidence that 2D HI is important for the shock wave formation { for a given A , an ordered shock is only formed when W exceeds a critical value. Second, it shows that the density shock formation is robust over a large plateau of W and A . Note that the innitely large limit of W is singular. For suciently large A , W large but nite means that the sidewalls suppress density uctuation in the y-direction, ensuring an approximately homogeneous distribution of swimmers in that direction. Whereas for W =1, this constraint disappears and the swimmers are free to escape to innity in y-direction, causing large density uctuation in the y-direction and suppression of the density shock wave in the x-direction. We use the shock order parameter to capture the the transition from subsonic to supersonic density shock wave via increasing the strength of external background ow, see Fig. 7.8. The sharp discontinuity in the shock order parameter at V = c 114 1 1.5 2 2.5 3 −1 −0.5 0 0.5 1 V < S > V > c V < c t → ∞ V = c Figure 7.8: Mean value of the shock order parameter showing the transition of density shock waves from subsonic to supersonic via increasing background ow velocity V, averaged over several sample simulations. Parameters used: W = 60, = 0:5, A = 0:3, = 1. indicates the transition of the density shock waves from subsonic to supersonic as the velocity of the background ow V increases. We also probe the long time dynamics of the shock wave in the periodic channel. The shock front persists for a long time, while with the rest of the swimmers ll the whole periodic channel, see Fig. 7.9. The shock front serves as a large localized density perturbation in a more or less homogenous medium and travels as a density band. Density band can also be observed when we consider an initial distribution of swimmers lling the whole channel with a large localized density perturbation of swimmers. In this case, the density perturbation travels like a band and persist for a long time, see Fig. 7.10. 115 0 200 400 0 0.4 0 0.4 0 200 400 200 400 0 0.4 0 ρ x x y 0 100 200 300 400 −20 0 20 0 100 200 300 400 −20 0 20 0 100 200 300 400 −20 0 20 y y ρ ρ t = 0 t = 1500 t = 5000 Figure 7.9: The shock front of microswimmers propagates in a periodic channel lled with swimmers eventually, and evolve as a traveling density band. Parameter values are: W = 40, L = 400, A =0.3, = 0:5, = 1, V = 3. 0 100 200 300 400 −20 0 20 0 100 200 300 400 −20 0 20 0 100 200 300 400 −20 0 20 0 200 400 0 0.6 0 200 0 0.6 400 0 200 400 0 0.6 y x x ρ y y ρ ρ t = 1000 t = 500 t = 0 Figure 7.10: Microswimmers homogeneously lling the ow channel: a large density perturbation propagates in the homogenous medium as a traveling density band. Parameter values are: W = 40,L = 400, A =0.3, = 0:5, = 1,V = 1:8, whereas A = 0:2 for the homogenous background density and A = 0:5 for the localized density perturbation. 7.5 Continuum equation model for density shock formation To further elucidate the eect of the complex interplay between the 2D HI, sidewall connement, and background ow on the density shock formation, we develop a 116 reduced continuum model that properly accounts for all these elements. Our model is distinct from the trac ow model used in [69] and Burgers' equation used in [8]. These equations oversimplify the eect of HI. They are only correct to leading order in the weak connement limit withaW and are not capable of explaining the transition of shock formation from subsonic to supersonic as the background ow velocity V increases. Here, we assume the swimmers' density is homogenous in they-direction and approximate the swimmers' number density as a 1D function c (x;t). The continuum number density c is related to via c ==(a 2 ). From (7.7), we develop a 1D continuum model (recall that p = sgn()) @ c @t + @ @x h (Up +V +u) c i = 0: (7.22) The 2D HI are captured by the term u(x;t) u(x;t) = 1 W2a Re ZZZ w c (z;z 0 ) c (x 0 ;t)dydy 0 dx 0 (7.23) where x 0 is integrated from 0 to L and y, y 0 fromW=2 +a to W=2a, while the complex conjugate dipolar ow eld w c is given by w c = d a 2W 2 n csch 2 h 2W (zz 0 ) i sech 2 h 2W (zz 0 ) io : (7.24) We account for the steric repulsion between the swimmers via the excluded area jxx 0 j 2 +jyy 0 j 2 < 4a 2 , inside which the dipolar contribution is removed. The 117 eect of the excluded area is also re ected in the upper and lower limits of the y-integral in (7.10), indicating that no swimmers can penetrate the sidewall bound- aries. We solve the continuum equation (7.22) numerically with periodic boundary conditions using nte dierence method, with forward Euler scheme in time and rst-order upwind scheme in space. Note that the density function c in (7.23) does not depend on y, which enables us to pre-calculate the y-integral and reduce the integration from a triple integral to a single integral. This property largely reduces the computational cost of solving (7.22). The numerical solutions, superimposed in red in Figure 7.4, show excellent agree- ment between the discrete particle simulations and the continuum model, including successfully capturing the transition in the location of the shock formation as the background ow increases and the time scale of wave dispersion. 7.6 Control density distribution and population speed This physical understanding of the emergent behavior and its dependence on the swimmers' motility properties and intensity of background ow can be exploited to devise novel mechanisms for controlling populations of microswimmers in externally- driven ow channels. Two species of microswimmers of dierent motility properties subject to the same background ow V get segregated in space in nite time, see Fig. 7.11. This is particularly useful for biomedical applications such as sorting of 118 cells in ow channels. Another compelling scenario is to control the throughput, i.e., the speed of a cell population through a ow channel, by careful manipulation of the background owV . One could also control the prole of the population density distribution. For example, a Gaussian density prole can be obtained and main- tained from an initially uniform distribution by properly tailoring a time-dependent background owV , see Fig. 7.12. To achieve this, we introduce a weak background ow with V = 1 initially, which corresponds to the case of d > 0, and a density shock is formed at the back of the microswimmers. Then at t = 1000, we switch the background ow to strong intensity, that is, we set V = 3 and thus d < 0. In this case, the swimmers at the back of the wave tend to push themselves towards the front of the wave, whereas the swimmers at the front of the wave are pulled to the back, and therefore results in a Gaussian density prole at t = 1500. The Gaussian density prole is then maintained by a background ow of fast switching intensity such that the polarization vanishes on average, that ish d i t = 0. The intensity of background ow is switched between two values, which are V = 1 and V = 3, and the period of switching is 100 dimensionless time units. The Gaussian density prole keeps qualitatively unchanged, as shown in the prole at t = 2000. 7.7 Conclusions This work provides the rst systematic study of the emergence of compression shock waves in populations of mircoswimmers driven in narrow ow channels. Our results are consistent with experimental observations on driven droplets and self-propelled 119 −40 0 40 −40 0 40 0 200 400 600 800 1000 −40 0 40 x y V = 1 μ = 0.5 υ = −1 μ = 0.5 υ = 1 μ = 0.25 υ = 1 μ = 0.5 υ = 1 t = 0 t = 150 t = 750 (a) (b) (c) y y Figure 7.11: Segregation of two species of microswimmers that are mixed homoge- neously, with initial conditions given in (a). The two species of swimmers are of distinct (b) translational and (c) rotational motility, respectively. Parameter values are: W = 80, L = 800, A = 0:2 for both species of swimmers. 0 0.2 0.4 −40 0 40 0 250 500 750 x ρ 1000 y t = 0 t = 1000 t = 1500 t = 2000 Figure 7.12: The strength of the background ow is controlled to obtain a Gaussian density distribution of microswimmers from an initially uniform density distribu- tion. Parameter Values are: W = 80,L = 800, A = 0:3, = 0:5, = 1. Snapshots of the swimmers at dierent time are shown on the same axis for illustration. colloids [7, 8, 11, 12], but we go beyond these observations to report an unobserved transition from `subsonic' to `supersonic' compression shocks as the intensity of the background ow increases. This transition from subsonic to supersonic shocks is only observed in self-propelled swimmers driven by a background ow. The phys- ical mechanisms underlying this transition are elucidated via discrete simulations 120 and a continuum model that properly captures the eects of HI and geometric connement. A rigorous analysis of the onset and propagation of instability in the continuum model is deferred to future work. These results, albeit in the context of a simplied model, have profound implications on developing a physics-based, quan- titative framework for the design and control of particle-laden ows in micro uidic channels. 121 Chapter 8 Conclusions The collective behavior of microswimmers is an active area of research. Most of these works have focused on the microswimming in an unbounded three-dimensional uid ow or in uid systems subject to weak connement. However, motivated by recent technological advances in producing and manipulating large ensembles of particles in micro uidic devices, there is growing interest in the problem of collective dynamics of swimming and motile particles strongly conned in quasi two-dimensional ge- ometries such as Hele-Shaw cells. Hele-Shaw connement changes drastically the hydrodynamics of microswimmer. A microswimmer conned in a Hele-Shaw cell has distinct hydrodynamic signature, mobility properties and orientation dynamics compared to its unconned counterpart. In this thesis, we propose a set of novel mathematical and computational mod- els to investigate the collective dynamics of conned microswimmers in various Hele-Shaw geometries. In particular, we focus on the transitions in the collective behavior due to the changes in the system parameters. Here, we summarize and 122 highlight several important results of this thesis, and provide the future outlook of the research. Role of agellar activity We extended the head-tail asymmetric swimmer model proposed by Brotto et al. [13] to take into account the eect of agellar activity. We showed that a conned mi- croswimmer with vigorously beating agella behaves like a large tail swimmer and tends to align itself with the local ow velocity, whereas a swimmer with weakly beating agella behaves like a large head swimmer and tends to reorient to the opposite of the local ow. It turns out that, this dierence in orientation dynamics due to variation in agellar activity is the hallmark of transitions in the collective behavior for a population of conned microswimmers. We show that decreasing ag- ellar activity induces a hydrodynamically triggered transition across three distinct modes of collective behavior, from the chaotic turbulentlike swimming to global ori- entational order to stable cluster. These results suggest that the interplay between agellar activity and hydrodynamic interactions provides a physical mechanism for coordinating collective behaviors in conned bacteria, with potentially profound implications on processes such as molecular diusion and transport of oxygen and nutrients that mediate transitions in the bacteria physiological state. Spontaneous global order in circular connement Inspired by recent experiments which demonstrated the self-organization of a pop- ulation of bacteria within a cylindrical droplet into a coherent vortex structure, 123 we considered the interplay between agellar activity, Hele-Shaw and circular con- nement on the emergent global patterns of the microswimmers. We proposed a model that properly captures the hydrodynamic eect of the radial wall and showed that decreasing agellar activity induces a hydrodynamically triggered transition in conned microswimmers from swirling to global vortex circulation to boundary ag- gregation and clustering. These results highlight that the complex interplay between connement, agellar activity, and hydrodynamic ows in concentrated suspensions of microswimmers could lead to a plethora of global patterns that are dicult to predict from geometric consideration alone. Emergence of density shock waves Motile and driven particles conned in micro uidic channels exhibit interesting emergent behavior from propagating density bands to density shock waves. A deeper understanding of the physical mechanisms responsible for these emergent structures is relevant to a number of physical and biomedical applications. We studied the formation of density shock waves in the context of microswimmers con- ned in a narrow channel and subject to a uniform external ow. The density shock waves emerge as a result of the dipolar interactions between the swimmers. These density shock waves exhibit a transition from `subsonic' with compression at the back to `supersonic' with compression at the front of the population as the intensity of the external ow increases. Interestingly, the location of the density shock in the two-dimensional system is completely opposite to the one obtained from the one-dimensional system. The underlying physical mechanisms for the 124 shock formation is attributed to the two-dimensional hydrodynamic interactions and sidewalls connement, and is conrmed by a novel quasilinear wave model that properly captures the dependence of the shock formation on the external ow. This physical understanding of the emergent behavior and its dependence on the swim- mers' motility properties and intensity of background ow can be exploited to devise novel mechanisms for controlling populations of microswimmers in externally-driven ow channels. In particular, we present several useful examples for the control of microswimmers, including segregation of microswimmers of dierent mobility prop- erties and the tailoring of emergent density distribution of the microswimmers into a desired prole. These ndings have potentially profound implications on various processes in industry and biotechnology such as the transport and sorting of cells in ow channels. Future work The models developed in this thesis have successfully captured the transition of collective behavior of conned microswimmers in dierent Hele-Shaw geometries. They help us to obtain simple physical pictures that explain the emergence of these collective behaviors. However, several research questions are left to be answered, including the ow structure, diusion and transport of passive particles in the col- lective ow, and the eective viscosity of the uid. Most of the work in this thesis is performed using the discrete particle approach. In the future, we will extend our analysis using the continuum approach. In fact, there exists a mismatch between the collective modes of conned microswimmers 125 predicted from the linear instability of the continuum model [13] and the results obtained from discrete particle simulations ([69] and Chapter 5). This is not sur- prising, as the linear stability analysis obtained from the continuum model is valid for the dilute populations of particle, whereas the particle simulations are performed on the semi-dilute regime. Future work is devoted to revisit the continuum model to discern the nonlinear eects on the emergent collective dynamics. Fully nonlinear simulations on the continuum model that captures the eects of excluded area will be performed. Density waves such as sound waves and shock waves are observed in populations of conned microswimmers, as shown in Chapter 7. The propagation of the density waves is attributed to the dipolar interaction between the swimmers. Future exten- sions of this work will include the study of the propagation of sound waves in the form of `phonons' in active micro uidic crystals, that is, microswimmers positioned in regular lattices. We will focus on the wave-wave interactions and the instability phenomena due to the wave disturbances. 126 References [1] M. Abramowitz and I. A. Stegun. Handbook of mathematical functions: with formulas, graphs, and mathematical tables. Number 55. Courier Corporation, 1964. [2] I. S. Aranson, A. Sokolov, J. O. Kessler, and R. E. Goldstein. Model for dynamical coherence in thin lms of self-propelled microorganisms. Physical Review E, 75(4):040901, 2007. [3] A. Baskaran and M. C. Marchetti. Statistical mechanics and hydrodynamics of bacterial suspensions. Proceedings of the National Academy of Sciences, 106(37):15567{15572, 2009. [4] G. K. Batchelor. An Introduction to Fluid Dynamics. Cambridge University Press, 1967. [5] T. Beatus, R. Bar-Ziv, and T. Tlusty. Anomalous micro uidic phonons in- duced by the interplay of hydrodynamic screening and incompressibility. Phys- ical Review Letters, 99(12):124502, 2007. [6] T. Beatus, R. H. Bar-Ziv, and T. Tlusty. The physics of 2d micro uidic droplet ensembles. Physics Reports, 516(3):103{145, 2012. [7] T. Beatus, T. Tlusty, and R. Bar-Ziv. Phonons in a one-dimensional micro u- idic crystal. Nature Physics, 2(11):743{748, 2006. [8] T. Beatus, T. Tlusty, and R. Bar-Ziv. Burgers shock waves and sound in a 2d micro uidic droplets ensemble. Physical Review Letters, 103(11):114502, 2009. [9] A. P. Berke, L. Turner, H. C. Berg, and E. Lauga. Hydrodynamic at- traction of swimming microorganisms by surfaces. Physical Review Letters, 101(3):038102, 2008. 127 [10] S. Bhattacharya, J .B lawzdziewicz, and E. Wajnryb. Far-eld approximation for hydrodynamic interactions in parallel-wall geometry. Journal of Compu- tational Physics, 212(2):718{738, 2006. [11] A. Bricard, J.-B. Caussin, D. Das, C. Savoie, V. Chikkadi, K. Shitara, O. Chepizhko, F. Peruani, D. Saintillan, and D. Bartolo. Emergent vortices in populations of colloidal rollers. Nature Communications, 6:7470, 2015. [12] A. Bricard, J.-B. Caussin, N. Desreumaux, O. Dauchot, and D. Bartolo. Emer- gence of macroscopic directed motion in populations of motile colloids. Nature, 503(7474):95{98, 2013. [13] T. Brotto, J.-B. Caussin, E. Lauga, and D. Bartolo. Hydrodynamics of con- ned active uids. Physical Review Letters, 110(3):038101, 2013. [14] I. Buttinoni, J. Bialk e, F. K ummel, H. L owen, C. Bechinger, and T. Speck. Dynamical clustering and phase separation in suspensions of self-propelled colloidal particles. Physical Review Letters, 110(23):238301, 2013. [15] M. E. Cates. Diusive transport without detailed balance in motile bacteria: does microbiology need statistical physics? Reports on Progress in Physics, 75(4):042601, 2012. [16] J. B. Caussin, A. Solon, A. Peshkov, H. Chat e, T. Dauxois, J. Tailleur, V. Vitelli, and D. Bartolo. Emergent spatial structures in ocking models: a dynamical system insight. Physical Review Letters, 112(14):148102, 2014. [17] N. Champagne, E. Lauga, and D. Bartolo. Stability and non-linear response of 1d micro uidic-particle streams. Soft Matter, 7(23):11082{11085, 2011. [18] H. Chat e, F. Ginelli, G. Gr egoire, and F. Raynaud. Collective motion of self-propelled particles interacting without cohesion. Physical Review E, 77(4):046113, 2008. [19] X. Chen, X. Yang, M. Yang, and H. P. Zhang. Dynamic clustering in suspen- sion of motile bacteria. Europhysics Letters, 111(5):54002, 2015. [20] A. T. Chwang and T. Wu. Hydromechanics of low-reynolds-number ow. part 2. singularity method for stokes ows. Journal of Fluid Mechanics, 67(04):787{815, 1975. [21] L. H. Cisneros, R. Cortez, C. Dombrowski, R. E. Goldstein, and J. O. Kessler. Fluid dynamics of self-propelled microorganisms, from individuals to concen- trated populations. Experiments in Fluids, 43(5):737{753, 2007. 128 [22] L. H. Cisneros, J. O. Kessler, S. Ganguly, and R. E. Goldstein. Dynamics of swimming bacteria: Transition to directional order at high concentration. Physical Review E, 83(6):061907, 2011. [23] R. Coquereaux, A. Grossmann, and B. E. Lautrup. Iterative method for calculation of the weierstrass elliptic function. IMA Journal of Numerical Analysis, 10(1):119{128, 1990. [24] A. Costanzo, R. Di Leonardo, G. Ruocco, and L. Angelani. Transport of self-propelling bacteria in micro-channel ow. Journal of Physics: Condensed Matter, 24(6):065101, 2012. [25] I. D. Couzin, J. Krause, N. R. Franks, and S. A. Levin. Eective leadership and decision-making in animal groups on the move. Nature, 433(7025):513{ 516, 2005. [26] I. D. Couzin, J. Krause, R. James, G. D. Ruxton, and N. R. Franks . Collective memory and spatial sorting in animal groups. Journal of Theoretical Biology, 218(1):1{11, 2002. [27] B. Cui, H. Diamant, B. Lin, and S. Rice. Anomalous hydrodynamic in- teraction in a quasi-two-dimensional suspension. Physical Review Letters, 92(25):258301, 2004. [28] N. C. Darnton, L. Turner, S. Rojevsky, and H. C. Berg. Dynamics of bacterial swarming. Biophysical Journal, 98(10):2082{2090, 2010. [29] M. E. Davey and G. A. O'toole. Microbial biolms: from ecology to molecular genetics. Microbiology and Molecular Biology Reviews, 64(4):847{867, 2000. [30] R. D. Deegan, O. Bakajin, T. F. Dupont, G. Huber, S. R. Nagel, and T. A. Witten. Capillary ow as the cause of ring stains from dried liquid drops. Nature, 389(6653):827{829, 1997. [31] B. Delmotte, E. Keaveny, F. Plouraboue, and E. Climent. Large-scale simula- tion of steady and time-dependent active suspensions with the force-coupling method. Journal of Computational Physics, 302:524{547, 2015. [32] N. Desreumaux, J.-B. Caussin, R. Jeanneret, E. Lauga, and D. Bartolo. Hy- drodynamic uctuations in conned particle-laden uids. Physical Review Letters, 111(11):118301, 2013. [33] N. Desreumaux, N. Florent, E. Lauga, and D. Bartolo. Active and driven hydrodynamic crystals. The European Physical Journal E, 35(8):68, 2012. 129 [34] C. Dombrowski, L. Cisneros, S. Chatkaew, R. E. Goldstein, and J. O. Kessler. Self-concentration and large-scale coherence in bacterial dynamics. Physical Review Letters, 93(9):098103, 2004. [35] J. Dunkel, S. Heidenreich, M. B ar, and R. E. Goldstein. Minimal continuum theories of structure formation in dense active uids. New Journal of Physics, 15(4):045016, 2013. [36] J. Dunkel, S. Heidenreich, K. Drescher, H. H. Wensink, M. B ar, and R. E. Goldstein. Fluid dynamics of bacterial turbulence. Physical Review Letters, 110(22):228102, 2013. [37] A. C. Edmunds, L. F. Castiblanco, G. W. Sundin, and C. M. Waters. Cyclic di-gmp modulates the disease progression of erwinia amylovora. Journal of Bacteriology, 195(10):2155{2165, 2013. [38] A. A. Evans, T. Ishikawa, T. Yamaguchi, and E. Lauga. Orientational order in concentrated suspensions of spherical microswimmers. Physics of Fluids, 23:111702, 2011. [39] B. Ezhilan, A. A. Pahlavan, and D. Saintillan. Chaotic dynamics and oxygen transport in thin lms of aerotactic bacteria. Physics of Fluids, 24(9):091701, 2012. [40] B. Ezhilan and D. Saintillan. Transport of a dilute active suspension in pressure-driven channel ow. Journal of Fluid Mechanics, 777:482{522, 2015. [41] B. Ezhilan, M. J. Shelley, and D. Saintillan. Instabilities and nonlinear dy- namics of concentrated activesuspensions. Physics of Fluids, 25(7):070607, 2013. [42] N. Figueroa-Morales, G. Mi~ no, A. Rivera, R. Caballero, E. Cl ement, E. Alt- shuler, and A. Lindner. Living on the edge: transfer and trac of e. coli in a conned ow. Soft Matter, 11:6284{6293, 2015. [43] J. Gachelin, G. Mi~ no, H. Berthet, A. Lindner, A. Rousselet, and E. Cl ement. Non-newtonian viscosity of escherichia coli suspensions. Physical Review Let- ters, 110(26):268103, 2013. [44] J. P. Hernandez-Ortiz, C. G. Stoltz, and M. D. Graham. Transport and collective dynamics in suspensions of conned swimming particles. Physical Review Letters, 95(20):204501, 2005. 130 [45] J. P. Hernandez-Ortiz, P. T. Underhill, and M. D. Graham. Dynamics of conned suspensions of swimming particles. Journal of Physics: Condensed Matter, 21(20):204107, 2009. [46] J. Hill, O. Kalkanci, J. L. McMurry, and H. Koser. Hydrodynamic surface interactions enable escherichia coli to seek ecient routes to swim upstream. Physical Review Letters, 98(6):068101, 2007. [47] C. Hohenegger and M. J. Shelley. Stability of active suspensions. Physical Review E, 81(4):046311, 2010. [48] J.-P. Hulin, A.-M. Cazabat, E. Guyon, and F. Carmona. Hydrodynamics of Dispersed Media. Elsevier, New York, 1990. [49] T. Ishikawa, J. T. Locsei, and T. J. Pedley. Development of coherent struc- tures in concentrated suspensions of swimming model micro-organisms. Jour- nal of Fluid Mechanics, 615:401{431, 2008. [50] T. Ishikawa, J. T. Locsei, and T. J. Pedley. Fluid particle diusion in a semidi- lute suspension of model micro-organisms. Physical Review E, 82(2):021408, 2010. [51] T. Ishikawa and T. J. Pedley. Diusion of swimming model micro-organisms in a semi-dilute suspension. Journal of Fluid Mechanics, 588:437{462, 2007. [52] T. Ishikawa and T. J. Pedley. Coherent structures in monolayers of swimming particles. Physical Review Letters, 100(8):088103, 2008. [53] T. Ishikawa and T. J. Pedley. Dispersion of model microorganisms swimming in a nonuniform suspension. Physical Review E, 90(3):033008, 2014. [54] E. L. Jackson and H. Lu. Advances in micro uidic cell separation and ma- nipulation. Current Opinion in Chemical Engineering, 2(4):398{404, 2013. [55] G. B. Jeery. The motion of ellipsoidal particles immersed in a viscous uid. Proceedings of the Royal Society of London Series A, 102(715):161{179, 1922. [56] E. Kanso and A. C. H. Tsang. Dipole models of self-propelled bodies. Fluid Dynamics Research, 46:061407, 2014. [57] E. Kanso and A. C. H. Tsang. Pursuit and synchronization in hydrodynamic dipoles. Journal of Nonlinear Science, 25(5):1141{1152, 2015. 131 [58] F. C. Keber, E. Loiseau, T. Sanchez, S. J. DeCamp, L. Giomi, M. J. Bowick, M. C. Marchetti, Z. Dogic, and A. R. Bausch. Topology and dynamics of active nematic vesicles. Science, 345(6201):1135{1139, 2014. [59] M. J. Kim and K. S. Breuer. Enhanced diusion due to motile bacteria. Physics of Fluids, 16(9):L78{L81, 2004. [60] M. J. Kim and K. S. Breuer. Controlled mixing in micro uidic systems using bacterial chemotaxis. Analytical Chemistry, 79(3):955{959, 2007. [61] M. J. Kim and K. S. Breuer. Micro uidic pump powered by self-organizing bacteria. Small, 4(1):111{118, 2008. [62] D. L. Koch and G. Subramanian. Collective hydrodynamics of swimming microorganisms: Living uids. Annual Review of Fluid Mechanics, 43:637{ 659, 2011. [63] S. K ohler, V. Schaller, and A. R. Bausch. Collective dynamics of active cytoskeletal networks. PLoS ONE, 6(8):e23798, 08 2011. [64] S. K ohler, V. Schallerr, and A. R. Bausch. Structure formation in active networks. Nature Materials, 10(6):462{468, 2011. [65] H. Kurtuldu, J. S. Guasto, K. A. Johnson, and J. P. Gollub. Enhancement of biomixing by swimming algal cells in two-dimensional lms. Proceedings of the National Academy of Sciences, 108(26):10391{10395, 2011. [66] E. Lauga, W. R. DiLuzio, G. M. Whitesides, and H. A. Stone. Swimming in circles: motion of bacteria near solid boundaries. Biophysical Journal, 90(2):400{412, 2006. [67] E. Lauga and R. E. Goldstein. Dance of the microswimmers. Physics Today, 65(9):30, 2012. [68] E. Lauga and T. R. Powers. The hydrodynamics of swimming microorganisms. Reports on Progress in Physics, 72(9):096601, 2009. [69] A. Lefauve and D. Saintillan. Globally aligned states and hydrodynamic trac jams in conned suspensions of active asymmetric particles. Physical Review E, 89:021002, 2014. [70] R. Di Leonardo, L. Angelani, D. Dell'Arciprete, G. Ruocco, V. Iebba, S. Schippa, M. P. Conte, F. Mecarini, F. De Angelis, and E. Di Fabrizio. Bacterial ratchet motors. Proceedings of the National Academy of Sciences, 107(21):9541{9545, 2010. 132 [71] K. C. Leptos, J. S. Guasto, J. P. Gollub, A. I. Pesci, and R. E. Goldstein. Dynamics of enhanced tracer diusion in suspensions of swimming eukaryotic microorganisms. Physical Review Letters, 103(19):198103, 2009. [72] G. J. Li and A. M. Ardekani. Hydrodynamic interaction of microswimmers near a wall. Physical Review E, 90(1):013010, 2014. [73] N. Liron and S. Mochon. Stokes ow for a stokeslet between two parallel at plates. Journal of Engineering Mathematics, 10(4):287{303, 1976. [74] H. M. L opez, J. Gachelin, C. Douarche, H. Auradou, and E. Cl ement. Turning bacteria suspensions into super uids. Physical Review Letters, 115:028301, Jul 2015. [75] E. Lushi, R. E. Goldstein, and M. J. Shelley. Collective chemotactic dynamics in the presence of self-generated uid ows. Physical Review E, 86(4):040902, 2012. [76] E. Lushi and C. S. Peskin. Modeling and simulation of active suspensions containing large numbers of interacting micro-swimmers. Computers & Struc- tures, 122:239{248, 2013. [77] E. Lushi, H. Wioland, and R. E. Goldstein. Fluid ows created by swimming bacteria drive self-organization in conned suspensions. Proceedings of the National Academy of Sciences, 110(27):201405698, 2014. [78] M. C. Marchetti. Soft matter: Frictionless uids from bacterial teamwork. Nature, 525(7567):37{39, 2015. [79] M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao, and R. A. Simha. Hydrodynamics of soft active matter. Reviews of Modern Physics, 85(3):1143, 2013. [80] H. Masoud and M. J. Shelley. Collective surng of chemically active particles. Physical Review Letters, 112(12):128304, 2014. [81] N. H. Mendelson, A. Bourque, K. Wilkening, K. R. Anderson, and J. C. Watkins. Organized cell swimming motions in bacillus subtilis colonies: pat- terns of short-lived whirls and jets. Journal of Bacteriology, 181(2):600{609, 1999. [82] M. B. Miller and B. L. Bassler. Quorum sensing in bacteria. Annual Reviews in Microbiology, 55(1):165{199, 2001. 133 [83] S. Mishra, A. Baskaran, and M. C. Marchetti. Fluctuations and pattern formation in self-propelled particles. Physical Review E, 81(6):061916, 2010. [84] K. J. Morton, K. Loutherback, D. W. Inglis, O. K. Tsui, J. C. Sturm, S. Y. Chou, and R. H. Austin. Crossing micro uidic streamlines to lyse, label and wash cells. Lab on a Chip, 8(9):1448{1453, 2008. [85] K. J. Morton, K. Loutherback, D. W. Inglis, O. K. Tsui, J. C. Sturm, S. Y. Chou, and R. H. Austin. Hydrodynamic metamaterials: Microfabricated ar- rays to steer, refract, and focus streams of biomaterials. Proceedings of the National Academy of Sciences, 105(21):7434{7438, 2008. [86] M. A. Mulvey, J. D. Schilling, J. J. Martinez, and S. J. Hultgren. Bad bugs and beleaguered bladders: interplay between uropathogenic escherichia coli and innate host defenses. Proceedings of the National Academy of Sciences, 97(16):8829{8835, 2000. [87] F. J. Ndlec, T. Surrey, A. C. Maggs, and S. Leibler. Self-organization of microtubules and motors. Nature, 389(6648):305{308, 1997. [88] P. K. Newton. The dipole dynamical system. Discrete and Continuous Dy- namical Systems (Suppl.), 2005:692{699, 2005. [89] K. A. O'Neil. On the hamiltonian-dynamics of vortex lattices. Journal of Mathematical Physics, 30(6):1373{1379, June 1989. [90] A. A. Pahlavan and D. Saintillan. Instability regimes in owing suspensions of swimming micro-organisms. Physics of Fluids, 23(1):011901, 2011. [91] O. S. Pak and E. Lauga. Theoretical models in low-Reynolds number locomo- tion. Chapter in Low-Reynolds-Number Flows: Fluid-Structure Interactions C. Duprat and H. A. Stone (Eds.), Royal Society of Chemistry Soft Matter Series, 2015. [92] J. Palacci, S. Sacanna, A. P. Steinberg, D. J. Pine, and P. M. Chaikin. Living crystals of light-activated colloidal surfers. Science, 339(6122):936{940, 2013. [93] A. P. Petro, X.-L. Wu, and A. Libchaber. Fast-moving bacteria self-organize into active two-dimensional crystals of rotating cells. Physical Review Letters, 114(15):158102, 2015. [94] A. S. Popel and P. C. Johnson. Microcirculation and hemorheology. Annual Review of Fluid Mechanics, 37:43, 2005. 134 [95] B. Qian, H. Jiang, D. A. Gagnon, K. S. Breuer, and T. R. Powers. Minimal model for synchronization induced by hydrodynamic interactions. Physical Review E, 80(6):61919, 2009. [96] S. Rafa , L. Jibuti, and P. Peyla. Eective viscosity of microswimmer suspen- sions. Physical Review Letters, 104(9):098102, 2010. [97] S. Ramaswamy. The mechanics and statistics of active matter. Annual Review of Condensed Matter Physics, 1:323{345, 2010. [98] W. Richtering. Smart Colloidal Materials, volume 133. Springer, 2006. [99] J. A. Riell and R. K. Zimmer. Sex and ow: the consequences of uid shear for sperm{egg interactions. The Journal of Experimental Biology, 210(20):3644{3660, 2007. [100] S. Sabrina, M. Spellings, S. C. Glotzer, and K. J. M. Bishop. Coarsening dynamics of binary liquids with active rotation. Soft matter, 11(43):8409{ 8416, 2015. [101] D. Saintillan and M. J. Shelley. Orientational order and instabilities in sus- pensions of self-locomoting rods. Physical Review Letters, 99(5):058102, 2007. [102] D. Saintillan and M. J. Shelley. Instabilities and pattern formation in ac- tive particle suspensions: kinetic theory and continuum simulations. Physical Review Letters, 100(17):178103, 2008. [103] D. Saintillan and M. J. Shelley. Instabilities, pattern formation, and mixing in active suspensions. Physics of Fluids, 20(12):123304, 2008. [104] D. Saintillan and M. J. Shelley. Emergence of coherent structures and large- scale ows in motile suspensions. Journal of The Royal Society Interface, 9(68):571{585, 2012. [105] D. Saintillan and M. J. Shelley. Active suspensions and their nonlinear models. Comptes Rendus Physique, 14(6):497{517, 2013. [106] D. Saintillan and M. J. Shelley. Theory of active suspensions. Springer, 2015. [107] T. Sanchez, D. T. N. Chen, S. J. DeCamp, M. Heymann, and Z. Dogic. Spontaneous motion in hierarchically assembled active matter. Nature, 491(7424):431{434, 2012. 135 [108] V. Schaller, C. Weber, C. Semmrich, E. Frey, and A. R. Bausch. Polar patterns of driven laments. Nature, 467(7311):73{77, 2010. [109] B. Scharf. Real-time imaging of uorescent agellar laments of rhizobium lupini h13-3: agellar rotation and ph-induced polymorphic transitions. Jour- nal of Bacteriology, 184(21):5979{5986, 2002. [110] J. Schwarz-Linek, C. Valeriani, A. Cacciuto, M. E. Cates, D. Marenduzzo, A. N. Morozov, and W. C. K. Poon. Phase separation and rotor self-assembly in active particle suspensions. Proceedings of the National Academy of Sci- ences, 109(11):4052{4057, 2012. [111] I. Shani, T. Beatus, R. H. Bar-Ziv, and T. Tlusty. Long-range orientational order in two-dimensional micro uidic dipoles. Nature Physics, 10(2):140{144, 2014. [112] R. A. Simha and S. Ramaswamy. Hydrodynamic uctuations and instabilities in ordered suspensions of self-propelled particles. Physical Review Letters, 89(5):058101, 2002. [113] A. Sokolov and I. S. Aranson. Reduction of viscosity in suspension of swim- ming bacteria. Physical Review Letters, 103(14):148101, 2009. [114] A. Sokolov, I. S. Aranson, J. O. Kessler, and R. E. Goldstein. Concentration dependence of the collective dynamics of swimming bacteria. Physical Review Letters, 98(15):158102, 2007. [115] A. Sokolov, R. E. Goldstein, F. I. Feldchtein, and I. S. Aranson. Enhanced mixing and spatial instability in concentrated bacterial suspensions. Physical Review E, 80(3):031903, 2009. [116] S. E. Spagnolie and E. Lauga. Hydrodynamics of self-propulsion near a bound- ary: predictions and accuracy of far-eld approximations. Journal of Fluid Mechanics, 700:105{147, 2012. [117] M. Spellings, M. Engel, D. Klotsa, S. Sabrina, A. M. Drews, N. H. P. Nguyen, K. J. M. Bishop, and S. C. Glotzer. Shape control and compartmentalization in active colloidal cells. Proceedings of the National Academy of Sciences, 112(34):E4642{E4650, 2015. [118] M. A. Stremler. On relative equilibria and integrable dynamics of point vor- tices in periodic domains. Theoretical and Computational Fluid Dynamics, 24(1):25{37, 2010. 136 [119] M. A. Stremler and H. Aref. Motion of three point vortices in a periodic parallelogram. Journal of Fluid Mechanics, 392:101{128, 1999. [120] G. Subramanian and D. L. Koch. Critical bacterial concentration for the onset of collective swimming. Journal of Fluid Mechanics, 632:359{400, 2009. [121] Y. Sumino, K. H. Nagai, Y. Shitaka, D. Tanaka, K. Yoshikawa, H. Chat e, and K. Oiwa. Large-scale vortex lattice emerging from collectively moving microtubules. Nature, 483(7390):448{452, 2012. [122] A. Sun and J. Lahann. Dynamically switchable biointerfaces. Soft Matter, 5(8):1555{1561, 2009. [123] T. Surrey, F. N ed elec, S. Leibler, and E. Karsenti. Physical proper- ties determining self-organization of motors and microtubules. Science, 292(5519):1167{1171, 2001. [124] Georey Taylor. Analysis of the swimming of microscopic organisms. Pro- ceedings of the Royal Society of London Series A, 209(1099):447{461, 1951. [125] I. Theurkau, C. Cottin-Bizonne, J. Palacci, C. Ybert, and L. Bocquet. Dy- namic clustering in active colloidal suspensions with chemical signaling. Phys- ical Review Letters, 108(26):268303, 2012. [126] S. Thutupalli, R. Seemann, and S. Herminghaus. Swarming behavior of simple model squirmers. New Journal of Physics, 13(7):073021, 2011. [127] V. K. Tkachenko. On vortex lattices. Soviet Journal of Experimental and Theoretical Physics, 22:1282, 1966. [128] J. Toner, Y. Tu, and S. Ramaswamy. Hydrodynamics and phases of ocks. Annals of Physics, 318(1):170{244, 2005. [129] A. C. H. Tsang and E. Kanso. Dipole interactions in doubly periodic domains. Journal of Nonlinear Science, 23:971{991, 2013. [130] A. C. H. Tsang and E. Kanso. Flagella-induced transitions in the collective behavior of conned microswimmers. Physical Review E, 90:021001, 2014. [131] A. C. H. Tsang and E. Kanso. Circularly conned microswimmers exhibit multiple global patterns. Physical Review E, 91:043008, 2015. [132] A. C. H. Tsang and E. Kanso. Density shock waves in conned microswim- mers. Physical Review Letters, 116:048101, 2016. 137 [133] W. E. Uspal, H. B. Eral, and P. S. Doyle. Engineering particle trajectories in micro uidic ows using particle shape. Nature Communications, 4:2666, 2013. [134] T. Vicsek, A. Czir ok, E. Ben-Jacob, I. Cohen, and O. Shochet. Novel type of phase transition in a system of self-driven particles. Physical Review Letters, 75(6):1226, 1995. [135] T. Vicsek and A. Zafeiris. Collective motion. Physics Reports, 517(3):71{140, 2012. [136] I. D. Vladescu, E. J. Marsden, J. Schwarz-Linek, V. A. Martinez, J. Arlt, A. N. Morozov, D. Marenduzzo, M. E. Cates, and W. C. K. Poon. Filling an emulsion drop with motile bacteria. Physical Review Letters, 113(26):268101, 2014. [137] J. Vollmer, A. G. Vegh, C. Lange, and B. Eckhardt. Vortex formation by active agents as a model for Daphnia swarming. Physical Review E, 73:061924, 2006. [138] A. Walther and A. E. M uller. Janus particles: synthesis, self-assembly, phys- ical properties, and applications. Chemical Reviews, 113(7):5194{5261, 2013. [139] H. H. Wensink, J. Dunkel, S. Heidenreich, K. Drescher, R. E. Goldstein, H. L owen, and J. M. Yeomans. Meso-scale turbulence in living uids. Pro- ceedings of the National Academy of Sciences, 109(36):14308{14313, 2012. [140] H. H. Wensink, V. Kantsler, R. E. Goldstein, and J. Dunkel. Controlling active self-assembly through broken particle-shape symmetry. Physical Review E, 89(1):010302, 2014. [141] E. T. Whittaker and G. N. Watson. A course of modern analysis. Cambridge university press, 1996. [142] H. Wioland, F. G. Woodhouse, J. Dunkel, J. O. Kessler, and R. E. Goldstein. Connement stabilizes a bacterial suspension into a spiral vortex. Physical Review Letters, 110(26):268102, 2013. [143] F. G. Woodhouse and R. E. Goldstein. Spontaneous circulation of conned active suspensions. Physical Review Letters, 109(16):168105, 2012. [144] F. G. Woodhouse and R. E. Goldstein. Cytoplasmic streaming in plant cells emerges naturally by microlament self-organization. Proceedings of the Na- tional Academy of Sciences, 110(35):14132{14137, 2013. 138 [145] D. M. Woolley, R. F. Crockett, W. D. Groom, and S. G. Revell. A study of synchronisation between the agella of bull spermatozoa, with related obser- vations. Journal of Experimental Biology, 212(14):2215{2223, 2009. [146] X.-L. Wu and A. Libchaber. Particle diusion in a quasi-two-dimensional bacterial bath. Physical Review Letters, 84(13):3017, 2000. [147] K. Yeo, E. Lushi, and P. M. Vlahovska. Collective dynamics in a binary mixture of hydrodynamically coupled microrotors. Physical Review Letters, 114(18):188301, 2015. [148] H. P. Zhang, A. Be'er, E.-L. Florin, and H. L. Swinney. Collective motion and density uctuations in bacterial colonies. Proceedings of the National Academy of Sciences, 107(31):13626{13630, 2010. [149] H. P. Zhang, A. Be'er, R. S. Smith, E.-L. Florin, and H. L. Swinney. Swarming dynamics in bacterial colonies. Europhysics Letters, 87(4):48011, 2009. [150] A. Z ottl and H. Stark. Nonlinear dynamics of a microswimmer in poiseuille ow. Physical Review Letters, 108(21):218104, 2012. [151] A. Z ottl and H. Stark. Hydrodynamics determines collective motion and phase behavior of active colloids in quasi-two-dimensional connement. Physical Review Letters, 112(11):118101, 2014. 139 Appendices 140 Appendix A Derivation of ow velocity eld of potential dipole in a doubly periodic domain In this section, we derive the exact analytical closed-formed expression for a two dimensional potential dipole in a doubly periodic domain. The expression for a potential dipole in a doubly periodic domain is introduced in Chapter 3, and the interactions among a pair and a population of such potential dipoles are studied in detail in Chapters 4 and 5. A potential dipole can be constructed using fundamental singularity solutions of the Laplace's equation [56, 88]. We start by considering the derivation of a potential dipole in an unbounded domain. Consider two point vortices of equal and opposite strengths () placed a distance apart. The vortex of strength is referred to as the left vortex and its position is denoted byz ` whereas the position of the right 141 vortex is denoted by z r . In the complex plane (z = x + iy and i = p 1), the locations of the point vortices are given by z ` =z o + ie io 2 ; z r =z o ie io 2 ; (A.1) wherez o is the center of the two point vortices. The complex conjugate ow velocity eld induced by the point vortices has the form w(z) = i 1 zz r 1 zz ` : (A.2) Substitute (A.1) into (A.2), we obtain w(z) = e io (zz o + ie io )(zz o ie io ) : (A.3) Let the distance between the two point vortices be small, say, ! 0, but assume that the dipole strength is of order 1. Then from (A.3), we obtain the ow velocity eld of a potential dipole in the unbounded domain, w(z) = e io (zz o ) 2 ; (A.4) where the location of the dipole is denoted by z o , and o is the orientation of the dipole. 142 Now, we derive the expression of potential dipole in a doubly periodic domain. We follow similar analysis performed in unbounded domain and consider the so- lution of a pair of point vortices in a doubly periodic domain. The formulation of point vortex in a doubly periodic domain has been derived in several dierent studies [89, 118, 119, 127, 129]. Building upon these results, the ow velocity eld induced by a pair of doubly periodic point vortices are given by w(z) = i [ (zz r ;! 1 ;! 2 ) (zz ` ;! 1 ;! 2 ) + ! 1 ! 1 1 ! 1 (z r z ` ) (z r z ` ) ]: (A.5) where (z;! 1 ;! 2 ) is the Weierstrass -function, (z;! 1 ;! 2 ) = 1 z + X k;l 1 z kl + 1 kl + z 2 kl : kl = 2k! 1 + 2l! 2 ; k;l2Zf0g: (A.6) Here, k and l are signed integers. In (A.5), ! 1 and ! 2 are the half-periods of the doubly-periodic domain, 1 is the value of the Weierstrass -function at the half- period! 1 , and is the area of the rectangular domain and is given by the Legendre's relation, = 2i(! 1 ! 2 ! 1 ! 2 ). The Weierstrass -function, albeit having innite 143 simple poles located at z = 0 and z = kl , is not a periodic function but a quasi- periodic one. The third and forth terms in (A.5) are the correction terms that enforce the periodicity, that is, to satisfy the condition, w(z + kl ) =w(z); (A.7) with the quasi-periodic property of Weierstrass -function (z + kl ) =(z) + 2k(! 1 ) + 2l(! 2 ): (A.8) Now, we consider the small distance limit of the periodic vortex pairs to obtain the potential dipole solution, that is, set ! 0 and keeps =. It turns out that in (A.5), the term corresponding to the dierence between the Weierstrass-functions reduces to the Weierstrass elliptic function (z), and the last two terms in (A.5) vanish in the limiting process. The ow velocity eld for a potential dipole in a doubly periodic domain is therefore given by w(z) =(zz o ;! 1 ;! 2 )e io ; (A.9) where Weierstrass elliptic function (z) is (z;! 1 ;! 2 ) = 1 z + X k;l 1 (z kl ) 2 1 2 kl : kl = 2k! 1 + 2l! 2 ; k;l2Zf0g: (A.10) 144 Appendix B Evaluation of Weierstrass type elliptic functions In this section, we discuss the ecient numerical methods of evaluating the Weier- strass elliptic function (z). The Weierstrass elliptic function (z;! 1 ;! 2 ) takes the form, (z;! 1 ;! 2 ) = 1 z 2 + X k;l 1 (z kl ) 2 1 2 pq : kl = 2k! 1 + 2l! 2 ; k;l2Zf0g: (B.1) The parameters k and l are signed integers. ! 1 and ! 2 are the half-periods of the doubly periodic domain. Direction computation of double sum in (B.1) is slow in convergence. Instead of computating the values of Weierstrass elliptic function from (B.1) directly, we can compute its functional values based on a fast convergent innite series expansion. 145 The Weierstrass -function can be expressed by the \q-series" expansion [141], (z;! 1 ;! 2 ) = 1 ! 1 + 2! 1 2 csc 2 z 2! 1 2 2 ! 2 1 1 X k=1 kq 2k 1q 2k cos kz ! 1 ; q = exp i! 2 ! 1 ; (B.2) where 1 =(! 1 ) is the value of Weierstrass -function at the half period ! 1 . The half period value 1 can be calculated by 1 = 2 12! 1 2 2 ! 1 1 X k=1 kq 2k 1q 2k : (B.3) Note that 1 depends only on the cell size of the periodic domain, thus only needs to be evaluated once in a numerical simulation. The single innite sum in (B.2) is much faster in convergence compared to the innite double sum in (B.1). The q-series expansion can also be applied to other Weierstrass type function, including the derivative of Weierstrass elliptic function 0 , Weierstrass -function and Weierstrass -function. The results are summarized as follow: 0 (z) = 3 4! 3 1 cot z 2! 1 csc 2 z 2! 1 + 2 3 ! 3 1 1 X k=1 k 2 q 2k 1q 2k sin kz ! 1 : (B.4) (z) = 1 z ! 1 + 2! 1 cot z 2! 1 + 2 ! 1 1 X k=1 q 2k 1q 2k sin kz ! 1 : (B.5) 146 log((z)) = log 2! 1 + 1 z 2 2! 1 + log sin z 2! 1 + 4 1 X k=1 q 2k k(1q 2k ) sin 2 kz 2! 1 : (B.6) The Weierstrass-function and the Weierstrass-function are related to the Weier- strass -function through (z) = 0 (z); d dz log((z)) =(z): (B.7) Note that the summation terms in (B.2) and (B.4) do not converge in large k, therefore proper truncations of the summation terms is necessary. For numerical implementation, the innite sums in (B.2), (B.4) { (B.6) are truncated when the relative change in the values of the summations are smaller than the user-dened threshold. In addition to computing Weierstrass -function using q-series expansion, an even faster iterative method can be developed based on the Laurent series and duplication formula of [23]. These formula can be found in a standard handbook of mathematical functions, such as Ref. [1]. Consider a pointz =z o near the origin, the rst three terms of the Laurent series for at z 0 are given by (z o ) = 1 z o 2 + g 2 20 z 2 o + g 3 28 z 4 o ; (B.8) 147 where g 2 and g 3 are the Weierstrass invariants, which can be calculated by g 2 =4(e 1 e 2 +e 1 e 3 +e 2 e 3 ) = 2(e 2 1 +e 2 2 +e 2 3 ); g 3 = 4e 1 e 2 e 3 = 4 3 (e 3 1 +e 3 2 +e 3 3 ): (B.9) where e 1 =(! 1 ), e 2 =(! 2 ) and e 3 =(! 1 ! 2 ) are constants. The duplication formula of Weierstrass -function is given by (2z) =2(z) + [6 2 (z)g 2 =2] 2 4[4 3 (z)g 2 (z)g 3 ] : (B.10) Now, we outline the procedure of the iterative method for calculation of. Con- sider a non-zero complex numberz that lies in a primary periodicity cell containing the origin. To evaluate (z), let z o =z=2 N , where N is the iteration number. The initial value of the the iteration, 0 , is dened as 0 = (z o ) and (z o ) is obtained from (B.8). The value of in the n-th step of the iteration is calculated by a recursive formula that is developed based on (B.10), n =2 n1 + (6 2 n1 g 2 =2) 2 4(4 3 n1 g 2 n1 g 3 ) : (B.11) The recursive formula (B.11) is applied repeatedly until n = N and N is a good approximation of the Weierstrass elliptic function (z). Numerical testings showed thatN = 6 can already acquire the maximum accuracy of double precision, that is, 148 relative error of iterationO(10 15 ). With (z) calculated, the derivatives of the Weierstrass -function can be computed using the following relations: [ 0 (z)] 2 = 4 3 (z)g 2 (z)g 3 ; 00 (z) = 6 02 (z)g 2 =2; 000 (z) = 12(z) 0 (z): (B.12) The iterative method can also be extended to calculate 0 , - and -function based on their Laurent series and duplication formula. For completeness, we present here the Laurent series and duplication formula for the Weierstrass type functions. For the Laurent series, (z) = 1 z 2 + 1 X k=2 c k z 2k2 ; 0 (z) = 2 z 3 + 1 X k=2 (2k 2)c k z 2k3 ; (z) = 1 z 1 X k=2 c k z 2k1 2k 1 ; (z) = 1 X m;n=0 a m;n g 2 2 m (2g 3 ) n z 4m+6n+1 (4m + 6n + 1)! ; (B.13) where coecients c k , a m;n are given by c k = 3 (2k + 1)(k 3) k2 X m=2 c m c km a m;n = 3(m + 1)a m+1;n1 + 16 3 (n + 1)a m2;n+1 1 3 (2m + 3n 1)(4m + 6n 1)a m1;n ; (B.14) 149 and c 2 =g 2 =20, c 3 =g 3 =28, a 0;0 = 0. The duplication formula for the Weierstrass type functions are given by (2z) =2(z) + 00 (z) 2 0 (z) 2 ; 0 (2z) = 4 4 (z) + 12(z) 02 (z) 00 (z) 003 (z) 4 03 (z) ; (2z) = 2(z) + 00 (z) 2 0 (z) ; (2z) = 0 (z) 4 (z); (3z) = 3(z) + 4 03 (z) 0 (z) 000 (z) 002 (z) ; (3z) = 02 (z) 9 (z)[(2z)(z)]; (B.15) where 0 (z), 00 (z), and 000 (z) are given by (B.12). 150
Abstract (if available)
Abstract
Microfluidics are revolutionizing modern research in disease diagnostics and drug development. For example, microfluidic devices for cell manipulation and sorting provide high-throughput methods for detection of infectious bacterial cells and offer great potential in cancer cell diagnostics. Great progress has been achieved in these directions, however, many questions concerning the design and control of these systems remain open. In particular, there is limited understanding of the collective behavior of motile cells in confined microfluidic environments, which makes the process of controlling these cells particularly challenging. In this thesis, we investigate the collective behavior of cell populations in the confined microfluidic geometries. ❧ We develop physics-based, mathematical and computational models of bacterial populations in confined microfluidic environment. We particularly focus on the hydrodynamic interactions among members of these bacterial populations. As such, our models account for the long-range hydrodynamic signature of bacterial cells under confinement and the short-range steric interactions among swimming bacteria and the microfluidic confinement. ❧ In one set of models, we investigate the effects of flagellar activity on the hydrodynamic interactions in bacterial populations in a doubly-periodic domain. We show that decreasing flagellar activity induces a hydrodynamically triggered transition in confined bacteria from swirling to aggregation. The swirling behavior is consistent with recent experimental findings, but our work goes beyond these findings to provide a predictive tool for systematically analyzing the collective bacterial behavior under various environmental conditions. In particular, our results show that varying the level of flagellar activity, which can be accomplished by varying environmental conditions, provide physical mechanisms for coordinating the collective behaviors in confined bacteria. This can in turn alter molecular microbial processes such as transport of oxygen and nutrients that mediate transitions in the bacteria physiological state. ❧ We then examine the effects of circular confinement on these transitions of collective behavior. We show as in the case of doubly-periodic domain that decreasing flagellar activity induces hydrodynamically triggered transitions, but from swirling to global circulation (vortex) to boundary aggregation and clustering. The spontaneous self-organization of bacterial populations into global vortex is reminiscent to recent experimental observations on circularly-confined populations of motile bacteria and rolling colloids. These results highlight that the geometry of the environments have significant influences on the emergent global patterns. ❧ Finally, we show that bacterial populations subject to increasing external flow in microfluidic channels exhibit an interesting transition of density shock waves from 'subsonic' to 'supersonic' modes, characterized by the shift in location of high population density from the back to the front of the group. These findings are consistent with experimental observations on driven droplets and self-propelled colloids. They also serve to guide the development of mechanisms for controlling the emergent density distribution and the average population speed of motile cells, with potentially profound implications on the transport and sorting of cells in microfluidic channels.
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
Synchronization in carpets of model cilia: numerical simulations and continuum theory
PDF
Passive flight in density-stratified fluid environments
PDF
Evaluating sensing and control in underwater animal behaviors
PDF
Microbial interaction networks: from single cells to collective behavior
PDF
Mechanical and flow interactions facilitate cooperative transport and collective locomotion in animal groups
PDF
Physics-based and data-driven models for bio-inspired flow sensing and motion planning
PDF
Multiscale modeling of cilia mechanics and function
PDF
Coordination in multi-muscle systems: from control theory to cortex
PDF
Passive rolling and flapping dynamics of a heaving Λ flyer
PDF
Traveling sea stars: hydrodynamic interactions and radially-symmetric motion strategies for biomimetic robot design
PDF
Mathematical characterizations of microbial communities: analysis and implications
PDF
A stochastic Markov chain model to describe cancer metastasis
PDF
Enabling human-building communication to promote pro-environmental behavior in office buildings
PDF
An experimental study of the elastic theory of granular flows
PDF
Understanding the formation and evolution of boundaries and interfaces in nanostructured metallic alloys
PDF
Modeling and simulation of circulating tumor cells in flow. Part I: Low-dimensional deformation models for circulating tumor cells in flow; Part II: Procoagulant circulating tumor cells in flow
PDF
A million-plus neuron model of the hippocampal dentate gyrus: role of topography, inhibitory interneurons, and excitatory associational circuitry in determining spatio-temporal dynamics of granul...
PDF
Speeding up trajectory planning for autonomous robots operating in complex environments
PDF
Quorum sensing interaction networks in bacterial communities
PDF
Computational tumor ecology: a model of tumor evolution, heterogeneity, and chemotherapeutic response
Asset Metadata
Creator
Tsang, Alan Cheng Hou
(author)
Core Title
Transitions in the collective behavior of microswimmers confined to microfluidic chips
School
Viterbi School of Engineering
Degree
Doctor of Philosophy
Degree Program
Mechanical Engineering
Publication Date
04/22/2016
Defense Date
03/07/2016
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
collective behavior,hydrodynamic interactions,microfluidics,microswimmers,OAI-PMH Harvest
Format
application/pdf
(imt)
Language
English
Contributor
Electronically uploaded by the author
(provenance)
Advisor
Kanso, Eva (
committee chair
), Haselwandter, Christoph (
committee member
), Kuhn, Peter (
committee member
), Luhar, Mitul (
committee member
), Newton, Paul (
committee member
)
Creator Email
alanchet@usc.edu,alantsangch@gmail.com
Permanent Link (DOI)
https://doi.org/10.25549/usctheses-c40-240423
Unique identifier
UC11278895
Identifier
etd-TsangAlanC-4352.pdf (filename),usctheses-c40-240423 (legacy record id)
Legacy Identifier
etd-TsangAlanC-4352.pdf
Dmrecord
240423
Document Type
Dissertation
Format
application/pdf (imt)
Rights
Tsang, Alan Cheng Hou
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Access Conditions
The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...
Repository Name
University of Southern California Digital Library
Repository Location
USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA
Tags
collective behavior
hydrodynamic interactions
microfluidics
microswimmers