Lecture 14 Notes: Low-Reynolds Number Limit

18.354J Nonlinear Dynamics II: Continuum Systems 14 Lecture 14 Low-Reynolds number limit In this section, we look at the limit of Re → 0 which is relevant to the construction of microfluidic devices and also governs the world of swimming microbes. Bacteria and eukaryotic cells achieve locomotion in a fluid through a self-induced change of shape: Escherichia coli propel themselves by rotating a helically shaped bundle of flagella, much like a corkscrew penetrating into a cork. Sperm cells move by inducing a wave-like deformation in a thin flagellum or cilium, whereas algae and other organisms move by beating two or more cilia in a synchronized manner (see slides). Because of their tiny size, these microswimmers operate at very low Reynolds number, i.e., inertial and turbulent effects are negligible15 . In this regime, swimming mechanisms are very different from employed by humans and other animals. In particular, any microbial swimming strategy must involve time-irreversible motion. Whilst moving through the liquid, a swimmer modifies the flow of the surrounding liquid. This can lead to an effective hydrodynamic interactions between nearby organisms, which can be attractive or repulsive depending on the details of the swimming mechanism. In reality, such deterministic forces are usually perturbed by a considerable amount of thermal or intrinsic noise, but we will neglect such Brownian motion effects here. 14.1 Stokes equations If the Reynolds number is very small, Re  1, the nonlinear NSEs (354) can be approximated by the linear Stokes equations 16 0 = ∇ · u, 2 0 = µ ∇ u − ∇p + f . (359a) (359b) The four equations (359) determine the four unknown functions (u, p). However, to uniquely identify such solutions, these equations must still be endowed with appropriate initial and boundary conditions, such as for example ( u(t, x) = 0, as |x| → ∞. (360) p(t, x) = p∞ , Note that, by neglecting the explicit time-dependent inertial terms in NSEs, the timedependence of the flow is determined exclusively and instantaneously by the motion of the boundaries and/or time-dependent forces as generated by the swimming objects. 14.2 Oseen’s solution Consider the Stokes equations (359) for a point-force f (x) = F δ(x). (361) 15 This is equivalent to larger animals swimming through a bath of treacle. More precisely, by replacing Eq. (354) with Eq. (359), it is assumed that for small Reynolds numbers ˜ Re(t, x) := |%(u · ∇)u|/(µ∇2 u) ' U L(%/µ)  1 one can approximate 16 % [∂t u + (u · ∇)u] − µ∇2 u ' −µ∇2 u The consistency of this approximation can be checked a posteriori by inserting the solution for u into the lhs. of Eq. (354) . In this case, the solution with standard boundary conditions (360) reads17 ui (x) = Gij (x) Fj , p(x) = Fj xj + p∞ , 4π|x|3 where the Greens function Gij is given by the Oseen tensor   xi xj 1 Gij (x) = δij + , 8πµ |x| |x|2 (362a) (362b) which has the inverse   xj xk G−1 (x) | | − , = 8πµ x δ jk jk 2|x|2 (363) as can be seen from Gij G−1 jk 14.3    xi xj xj xk = δij + δjk − |x|2 2|x|2 xi xj xj xk xi xk xi xk − = δik − + |x|2 2|x|2 |x|2 2|x|2 xi xk xi xk = δik − + 2|x|2 2|x|2 = δik . (364) Stokes’s solution (1851) Consider a sphere of radius a, which at time t is located at the origin, X(t) = 0, and moves at velocity U (t). The corresponding solution of the Stokes equation with standard boundary conditions (360) reads18      xj xi xj xi 3 a 1 a3 ui (t, x) = Uj δji + + δji − 3 2 , (365a) 4 |x| |x|2 4 |x|3 |x| 3 Uj x j p(t, x) = µa + p∞ . (365b) 2 |x|3 If the particle is located at X(t), one has to replace xi by xi −Xi (t) on the rhs. of Eqs. (365). Parameterizing the surface of the sphere by a = a sin θ cos φ ex + a sin θ sin φ ey + a cos θ ez = ai ei where θ ∈ [0, π], φ ∈ [0, 2π), one finds that on this boundary u(t, a(θ, φ)) = U , 3µ p(t, a(θ, φ)) = Uj aj (θ, φ) + p∞ , 2 a2 17 18 Proof by insertion. Proof by insertion. (366a) (366b) corresponding to a no-slip boundary condition on the sphere’s surface. The O(a/|x|)contribution in (365a) coincides with the Oseen result (362), if we identify F = 6π µa U . (367) The prefactor γ = 6π µa is the well-known Stokes drag coefficient for a sphere. The O[(a/|x|)3 ]-part in (365a) corresponds to the finite-size correction, and defining the Stokes tensor by   xj xi 1 a2 δji − 3 2 , (368) Sij = Gij + 24πµ |x|3 |x| we may rewrite (365a) as19 ui (t, x) = Sij Fj . 14.4 (369) Dimensionality We saw above that, in 3D, the fundamental solution to the Stokes equations for a point force at the origin is given by the Oseen solution ui (x) = Gij (x) Fj , p(x) = Fj xj + p∞ , 4π|x|3 (370a) where 1 Gij (x) = 8πµ |x|   xi xj δij + , |x|2 (370b) It is interesting to compare this result with corresponding 2D solution ui (x) = Jij (x)Fj , p= Fj xj + p∞ , 2π|x|2 x = (x, y) (371a) where     xi xj 1 |x| Jij (x) = −δij ln + 4πµ a |x|2 (371b) with a being an arbitrary constant fixed by some intermediate flow normalization condition. Note that (371) decays much more slowly than (370), implying that hydrodynamic interactions in 2D freestanding films are much stronger than in 3D bulk solutions. To verify that (371) is indeed a solution of the 2D Stokes equations, we first note that generally ∂j |x| = ∂j (xi xi )1/2 = xj (xi xi )−1/2 = xj |x| ∂j |x|−n = ∂j (xi xi )−n/2 = −nxj (xi xi )−(n+2)/2 = −n 19 For arbitrary sphere positions X(t), replace x → x − X(t). (372a) xj . |x|n+2 (372b) From this, we find Fj xj xi Fj Fi ∂i p = −2 = 2 4 2π|x| 2π|x| 2π|x|2   xj xi δij − 2 2 |x| (373) and ∂k Jij = = =     xi xj 1 |x| + ∂k −δij ln 4πµ a |x|2    xi xj 1 1 −δij ∂k |x| + ∂k 4πµ |x| |x|2    xi xj xk xj 1 xi xk . −δij 2 + δik 2 + δjk 2 − 2 4πµ |x| |x| |x| |x|4 To check the incompressibility condition, note that    xj xi xj xi 1 xi xi ∂i Jij = −δij 2 + δii 2 + δji 2 − 4πµ |x| |x| |x| 2|x|4   xj xj xj xj 1 = − 2 +2 2 + −2 2 2 4πµ |x| |x| |x| |x| = 0, (374) (375) which confirms that the solution (371) satisfies the incompressibility condition ∇ · u = 0. Moreover, we find for the Laplacian   xj xi xj xk ∂k xk xi ∂k ∂k Jij = −δij 2 + δik 2 + δjk 2 − 2 4πµ |x| |x| |x| |x|4          xj xi xj xk 1 xk xi = −δij ∂k + δik ∂k + δjk ∂k − 2∂k 4πµ |x|2 |x|2 |x|2 |x|4        δjk xj xk 1 δkk xk xk δik xi xk = − 2 −δij + δ − 2 + δ − 2 − ik jk 4πµ |x|2 |x|4 |x|2 |x|4 |x|2 |x|4   xi δjk xk xi xj δkk xi xj xk xk δik xj xk + + −4 2 |x|4 |x|4 |x|4 |x|6        δij xj xi δij xi xj 1 2 1 − = −δij − 2 + − 2 + − 2 4πµ |x|2 |x|2 |x|2 |x|4 |x|2 |x|4   xj xi xi xj xi xj xi xj 2 + +2 4 −4 4 |x|4 |x|4 |x| |x|   δij xj xi 1 (376) = −2 4 2πµ |x|2 |x| Hence, by comparing with (373), we see that indeed −∂i p + µ∂k ∂k ui = −∂i p + µ∂k ∂k Jij Fj = 0. (377) The difference between 3D and 2D hydrodynamics has been confirmed experimentally for Chlamydomonas algae. 14.5 Force dipoles In the absence of external forces, microswimmers must satisfy the force-free constraint. This simplest realization is a force-dipole flow, which provides a very good approximation for the mean flow field generated by an individual bacterium but not so much for a biflagellate alga. To construct a force dipole, consider two opposite point-forces F + = −F − = F ex located at positions x+ = ±`ex . Due to linearity of the Stokes equations the total flow at some point x is given by ui (x) = Γij (x − x+ ) Fj+ + Γij (x − x− ) Fj−   = Γij (x − x+ ) − Γij (x − x− ) Fj+ = [Γij (x − `ex ) − Γij (x + `ex )] Fj+ (378) where Γij = Jij in 2D and Γij = Gij in 3D. If |x|  `, we can Taylor expand Γij near ` = 0, and find to leading order    − ui (x) ' [Γij (x) − Γij (x)] − x+ Fj+ k ∂k Γij (x) − xk ∂k Γij (x) + = −2x+ k [∂k Γij (x)] Fj (379) 2D case Using our above result for ∂k Jij , and writing x+ = `n and F + = F n with |n| = 1, we find in 2D    xk+ xj xi xj xk xk xi ui (x) = − −δij 2 + δik 2 + δjk 2 − 2 Fj+ 2πµ |x| |x| |x| |x|4   xj nj nk xi xj xk nj F` xk n k xi −ni + ni = − + nk nk 2 − 2 |x|2 |x| |x|4 2πµ |x|2 and, hence, u(x) =  F`  2(n · x̂)2 − 1 x̂ 2πµ|x| (380) where x̂ = x/|x|. 3D case To compute the dipole flow field in 3D, we need to compute the partial derivatives of the Oseen tensor Gij (x) = 1 (1 + x̂i x̂j ) , 8πµ|x| x̂k = xk . |x| (381) Defining the orthogonal projector (Πik ) for x̂k by Πik := δik − x̂i x̂k , (382) we have xk = x̂k , |x| δik xk xi Πik = − = , 3 |x| |x| |x| 1 = − (x̂i Πnk + x̂k Πni ) , |x| ∂k |x| = ∂k x̂i ∂n Πik (383a) (383b) (383c) and from this we find ∂k Gij x̂k κ Gij + (Πik x̂j + Πjk x̂i ) |x| |x|2
κ −x̂ δ + x̂ δ + x̂ δ − 3x̂ x̂ x̂ . ij j i i j k ik jk k |x|2 = − = (384) Inserting this expression into (379), we obtain the far-field dipole flow in 3D u(x) =  F`  3(n · x̂)2 − 1 x̂. 2 4πµ|x| (385) Experiments show that Eq. (385) agrees well with the mean flow-field of a bacterium. Upon comparing Eqs. (380) and (385), it becomes evident that hydrodynamic interactions between bacteria in a free-standing 2D film are much longer-ranged than in a 3D bulk solution. This is a nice illustration of the fact that the number of available space dimensions can have profound effects on physical processes and interactions in biological systems. 14.6 Boundary effects The results in the previous section assumed an quasi-infinte fluid. Yet, many swimming cells and microorganisms operate in the vicinity of solid boundaries that can substantially affect the self-propulsion and the hydrodynamic interactions of the organisms. Before tackling finite boundaries geometries it is useful to recall how the terms in the HD equations can be rewritten in cylindrical coordinates. The full Navier-Stokes equations are written out completely in Acheson’s textbook for various coordinate systems. In the next part, we will summarize those terms that are most important for our further discussion of cylindrical confinements. 14.6.1 Reminder: Cartesian vs. cylindrical corodinates Cartesian coordinates In a global orthornormal Cartesian frame {ex , ey , ez }, the position vector is given by x = xex + yey + zez , and accordingly the flow field u(x) can be represented in the form u(x) = ux (x, y, z) ex + uy (x, y, z) ey + uz (x, y, z) ez . (386a) The gradient vector is given by ∇ = ex ∂x + ey ∂y + ez ∂z , (386b) and, using the orthonormality ej · ek = δjk , the Laplacian is obtained as ∆ = ∇ · ∇ = ∂x2 + ∂y2 + ∂z2 . (386c) One therefore finds for the vector-field divergence ∇ · u = ∂i ui = ∂x ux + ∂y uy + ∂z uz (386d) and the vector-Laplacian  ∂x2 ux + ∂y2 ux + ∂z2 ux ∆u = ∂i ∂i u =  ∂x2 uy + ∂y2 uy + ∂z2 uy  . ∂x2 uz + ∂y2 uz + ∂z2 uz  Cylindrical coordinates (386e) The local cylindrical basis vectors {er , eφ , ez } are defined by eφ = − sin φex + cos φey , er = cos φex + sin φey , φ ∈ [0, 2π) (387a) and they form a orthonormal system ej · ek = δjk , where now i, j = r, φ, z. The volume element is given by dV = r sin φ dr dφ dz. In terms of cylindrical basis system, the position vector x can be expressed as p x = rer + zez , r = x2 + y 2 (387b) (387c) and the flow field u(x) can be decomposed in the form u(x) = ur (r, φ, z) er + uφ (r, φ, z) eφ + uz (r, φ, z) ez . (387d) The gradient vector takes the form 1 ∇ = er ∂r + eφ ∂φ + ez ∂z , r (387e) yielding the divergence ∇·u = 1 1 ∂r (rur ) + ∂φ uφ + ∂z uz . r r (387f) The Laplacian of a scalar function f (r, φ, z) is given by ∇2 f = 1 1 ∂r (r∂r f ) + 2 ∂φ2 f + ∂z2 f r r (387g) and the Laplacian of a vector field u(r, φ, z) by ∇2 u = Lr er + Lφ eφ + Lz ez (387h) where Lr = Lφ = Lz = 1 1 2 1 ∂r (r∂r ur ) + 2 ∂φ2 ur + ∂z2 ur − 2 ∂φ uφ − 2 ur r r r r 1 1 2 2 1 2 ∂r (r∂r uφ ) + 2 ∂φ uφ + ∂z uφ + 2 ∂φ ur − 2 uφ r r r r 1 1 ∂r (r∂r uz ) + 2 ∂φ2 uz + ∂z2 uz r r (387i) (387j) (387k) Compared with the scalar Laplacian, the additional terms in the vector Laplacian arise from the coordinate dependence of the basis vectors. Similarly, one finds that, the r-component of (u·∇)u is not simply (u·∇)ur , but instead 1 er · [(u · ∇)u] = (u · ∇)ur − uφ2 . r (388) Physically, the term u2φ /r corresponds to the centrifugal force, and it arises because u = ur er + uφ eφ + uz ez and some of the unit vectors change with φ (e.g., ∂φ eφ = −er ). 14.6.2 Hagen-Poiseuille flow To illustrate the effects of no-slip boundaries on the fluid motion, let us consider pressure driven flow along a cylindrical pipe of radius R pointing along the z-axis. Assume that the flow is rotationally symmetric about the z-axis and constant in z-direction, u = uz (r)ez , p 2 2 where r = x + y is the distance from the center. For such a flow, the incompressibility condition ∇ · u = 0 is automatically satisfied, and the Stokes equation in cylindrical coordinates (r, φ, z) reduces to 0 = −∂z p + µ ∂r (r∂r uz ). r (389) Integrating twice over r, the general solution uz of this equation can be written as uz (r) = 1 (∂z p)r2 + c1 ln r + c2 , 4µ (390) where c1 and c2 are constants to be determined by the boundary conditions. For a no-slip boundary with uz (R) = 0 and finite flow speed at r = 0, one then finds uz (r) = − 1 (∂z p)(R2 − r2 ). 4µ (391) If we assume a linear pressure difference ∆P = P (L) − P (0) over a length L, then simply p(z ) = [P (L) − P (0)] z L ⇒ ∂z p = − P (0) − P (L) . L (392) The flow speed is maximal at center of the pipe u+ z = P (0) − P (L) 2 R 4µL (393) and the average transport velocity is uz = 1 πR2 Z R uz (r) 2πrdr = 0.5u+ z. (394) 0 Note that, for fixed pressure difference and channel length, the transport velocity uz decreases quadratically with the channel radius, signaling that the presence of boundaries can substantially suppress hydrodynamic flows. To illustrate this further, we next consider a useful approximation that can help to speed up numerical simulations through an effective reduction from 3D to 2D flow. 14.6.3 Hele-Shaw flow Consider two quasi-infinite parallel walls located at z = 0 and z = H. This setting is commonly encountered in experiments that study microbial swimming in flat microfluidic chambers. Looking for a 2D approximation of the Stokes equation, we may assume constant pressure along the z-direction, p = P (x, y), and neglect possible flow components in the vertical direction, uz = 0. Furthermore, using the above results for Hagen-Poiseuille flow as guidance, we can make the ansatz u(x, y, z) = 6z(H − z) 6z(H − z) U (x, y), [Ux (x, y)ex + Ux (x, y)ey ] ≡ 2 H H2 (395) corresponding to a parabolic flow profile in the vertical direction that accounts for no-slip boundaries at the walls; in particular, in the mid-plane 3 u(x, y, H/2) = U (x, y). 2 (396) We would like to obtain an effective equation for the effective 2D flow U (x, y). This can be achieved by inserting ansatz (395) into the Stokes equations and subsequently averaging along the z-direction20 , yielding 0 = −∇P + µ∇2 U − κU 0=∇·U , (397) where κ = 12µ/H 2 and ∇ is now the 2D gradient operator. Note that compared with unconfined 2D flow in a free film, the appearance of the κ-term leads to an exponential damping of hydrodynamic excitations. This is analogous to the exponential damping in the Yukawa-potential (mediated by massive bosons) compared to a Coloumb potential (mediated by massless photons). That ist, due to the presence of the no-slip boundaries, effective 2D hydrodynamic excitations acquire an effective mass ∝ 1/H 2 . 20 That is, by taking the integral (1/H) RH 0 dz of both sides.