2.4. Colloidal particles
2.4.6. Thermal Motion
We close this chapter by discussing how to include thermal noise in inertial microflu-idics. For most systems in inertial microfluidics, noise plays a secondary role compared to deterministic motion as we discuss below. However, it is still fruitful to include. First, stochastic fluctuations become important close to equilibrium positions where the de-terministic force vanishes. In particular, we include thermal noise in the description of feedback control in Sect. 5.2 to be able to observe the instability of the channel center.
Furthermore, even though thermal noise is a secondary effect it still is real in all physical systems. Finally, the method of multi-particle collisions dynamics, introduced in the next chapter, includes thermal noise. To interpret the simulation results, we have to understand thermal effects.
The forces the fluid exerts on colloidal particles originate from countless microscopic collisions between fluid molecules and the particles [83]. In addition to the forces de-scribed by the continuum equations, the collisions manifest themselves in an erratic motion of the colloid first described by Brown [83]. As we discuss below, the velocities in lateral and axial direction differ one to two orders in magnitude. We take this fact into account by including the thermal noise only in lateral direction. Therefore, we first con-sider exclusively the lateral coordinates x⊥ = (x, y) with gradient ∇⊥ = (∂/∂x, ∂/∂y) and comment on the axial coordinate z at the end of this section.
To include thermal noise, we add a stochastic forceη(t) to the equation of motion (2.27) for the particle. The corresponding equation of motion, known as Langevin equation, reads [83, 96]
M d2
dt2x⊥+ξd
dtx⊥=fext+η(t), (2.54) where we inserted the Stokes drag forceξv⊥ on the left hand side. We demand that the random noise has zero mean⟨ηi(t)⟩= 0 and is uncorrelated for different components and at different times,⟨ηi(t)ηj(t′)⟩= 2kBT ξδijδ(t−t′) [83]. The prefactor 2kBT ξ follows from the fluctuation dissipation theorem, discussed below. In a time interval ∆t, thermal noise changes the particle momentum by ∆p. The increments ∆p are independent Gaussian random variables with zero mean and variance 2kBT ξ∆t. We expect this form for many independent collisions with solvent molecules [83].
Similar to the reduction of the Navier-Stokes equations to the Stokes equation, the equation of motion for the colloid simplifies considerably in the friction-dominated regime.
If friction is large, the Langevin equation reduces to the overdamped Langevin equation [83, 96]
ξd
dtx⊥=fext+η(t). (2.55)
2. Basics of inertial microfluidics
To better understand the consequences of the stochastic force, we consider a free particle without external control forces (fext = 0) in the overdamped limit. In particular, we consider its displacement ∆x⊥(t) =x⊥(t)−x⊥(0) after time t. Since the stochastic force does not produce a mean drift, the ensemble average of the displacement vanishes [83, 96]
⟨∆x⊥(t)⟩= 0. (2.56)
On the other hand, the average of the squared displacement does not vanish. The correlations between displacements satisfy [83, 96]
⟨∆x⊥(t)⊗∆x⊥(t)⟩= 2Dt1, (2.57) where we introduced the diffusion constant D. The displacements in different directions are uncorrelated. Therefore, the trace of this expression is often considered, which results in⟨|∆x⊥(t)|2⟩= 2d Dt, wheredgives the number of lateral dimensions. For the Langevin equation, there exists a typical time scale M/ξ separating deterministic from diffusive motion [83]. For times t ≪ M/ξ, the particle moves with its initial velocity. For times t ≫M/ξ, the particle loses the memory of its initial velocity and diffuses [83].
The temperature and the translational friction coefficient determine the diffusion con-stant via the Einstein relation [83]
D= kBT
ξ . (2.58)
The Einstein relation is an example of the so-called fluctuation-dissipation theorem, which relates the equilibrium fluctuations (diffusion) to the response of the system to small forces (friction) [96, 97].
To estimate the relative importance of stochastic and deterministic forces, we intro-duce the P´eclet number as the ratio of potential energy to thermal energy. For inertial microfluidics, we estimate the potential energy by the product of inertial lift force flift and channel width 2w. Then, the P´eclet number takes the form
Pe = potential energy
thermal energy = 2wflift
kBT , (2.59)
In this thesis, the inertial lift force assumes typical valuesflift ≈1nN. For room temper-ature and a channel with width 2w = 20µm, the P´eclet number assumes typical values Pe ≈106. Hence, thermal noise plays only a role, when the inertial lift forces are small, for example, close to equilibrium positions.
In the analysis so far, we restricted ourselves to the first two moments of position in the case of free diffusion. To gain further insight, we have to consider the full probability
36
2.4. Colloidal particles
densityP(x⊥, t) to find a particle at positionx⊥ at timet. The probability density obeys the continuity equation [96, 97]
∂
∂tP(x⊥, t) +∇⊥·j(x⊥, t) = 0, (2.60) with probability flux j(x, t). The probability flux has two contributions: external forces result in deterministic motion, whereas random collisions result in diffusive motion. We write the total probability flux as [97]
j(x, t) =
ξ−1fext−D∇⊥
P(x⊥, t). (2.61)
The continuity equation with this form for the probability flux is also known as Smolu-chowski equation [96].
Alternatively, we interpret the divergence of the probability flux as the action of an operator L⊥ acting on the probability distribution [96]
∂
∂tP(x⊥, t) = L⊥P(x⊥, t), (2.62) where we defined the Smoluchowski operator
L⊥ =−∇⊥·
ξ−1fext−D∇⊥
. (2.63)
For conservative forces, we express the force in terms of a potential, fext = −∇⊥V. Then, the position in steady state follows the Boltzmann distribution [96]
P(x⊥) = 1 Z exp
− 1
kBTV(x⊥)
. (2.64)
Here, the partition function Z ensures that the distribution is normalized. A potential can always be constructed in one-dimensional systems or systems with radial symmetry.
When inertia is important, the Smoluchowski equation has to be extended to the Kramers equation. The steady-state distribution of the Kramers equation factorizes into position and velocity distributions. The position is distributed as before, whereas the velocity is distributed according to a Maxwell-Boltzmann distribution. The complete steady-state distribution of the Kramers equation reads [96]
P(x⊥,v⊥) = 1 Z exp
−1 2
m kBT|v⊥|2
exp
− 1
kBTV(x⊥)
. (2.65)
For bounded systems, we require boundary conditions for the Smoluchowski equation.
Here, we demand that the current through the boundary vanishes,nˆ·j = 0, wherenˆ is the normal to the boundary.
2. Basics of inertial microfluidics
Up to now, we discussed the Smoluchowski equation for a density evolving from an initial distribution. When we choose the initial distribution as a Dirac delta, the density becomes the probability P(x⊥,x′⊥|t) for the particle to move from its current position x′⊥to positionx⊥in time t. The transition probability satisfies a Smoluchowski equation too [96],
∂
∂tP(x,x′|t) =LP(x,x′|t), (2.66) where all spatial derivatives inLact on x. Since the particle was at position x′⊥att= 0 with certainty, the transition probability has to satisfy P(x⊥,x′⊥|t = 0) = δ(x⊥−x′⊥).
The no-flux boundary condition for the density translates into ˆ
n·
ξ−1fext−D∇⊥
P(x,x′|t) = 0 at the boundary. (2.67) We express the probability to find a particle at x⊥ at time t in terms of the transition probability and the distribution at t= 0 as
P(x⊥, t) =
dx′ P(x⊥,x′⊥|t)P(x′⊥,0). (2.68) We further use the transition probability to describe the evolution of expectation values in time. In particular, we consider the expectation value K(x⊥, t) of a function k(x⊥) at time t, when the particle starts at position x⊥(0). We express it as
K(x⊥, t) =
dx′⊥ P(x′⊥,x⊥|t)k(x′⊥), (2.69) where, we exchanged the primed and unprimed variable in the transition probability.
The time evolution of K(x, t) obeys the Kolmogorov backward equation [96, 98]
∂
∂tK(x⊥, t) =L+⊥K(x⊥, t), (2.70) where we introduced the adjoint operator
L+⊥ =
ξ−1fext+D∇⊥
·∇⊥. (2.71)
At the boundary the normal gradient of the expectation valueK has to vanish and we find [96]
ˆ
n·∇⊥K(x, t) = 0 at the boundary. (2.72) In the discussion up to this point, we only considered the lateral dimensions. We quantify the importance of thermal noise for the axial motion by the axial P´eclet number
38
2.4. Colloidal particles
Pez. Here, we use the drag force fdrag = ξv instead of the lift forces in Eq. (2.59). For typical velocities v ≈ 1m/s and typical particle radii a ≈ 5µm in inertial microfluidics, we find Pez ≈102Pe. Therefore, the axial motion is much less influenced by thermal noise and we assume the movement in axial direction to be deterministic. In microchannels, the particle translates with a velocityvz(x⊥) along the axial direction and the Langevin equation has to be supplemented by the equation
d
dtz =vz(x⊥). (2.73)
We include the axial motion in the Smoluchowski operator by adding a drift term along the axial direction. Then, the total Smoluchowski operator becomes
L=L⊥− ∂
∂zvz =−∇⊥·
ξ−1fext−D∇⊥
− ∂
∂zvz, (2.74) with its adjoint
L+ =L+⊥+vz ∂
∂z =
ξ−1fext+D∇⊥
·∇⊥+vz ∂
∂z. (2.75)
Mesoscopic simulations of fluid dynamics
Chapter 3
To correctly describe inertial migration, we have to solve the full Navier-Stokes equa-tions. Due their nonlinearity and the added complexity of moving boundaries from the colloidal particles, we choose to solve them numerically. In this chapter, we discuss the simulation methods used. First, we compare different methods to model fluid flow with a particular emphasis on mesoscopic simulation methods. Then, we describe the two particle-based mesoscopic methods used in this work: multi-particle collision dynamics and the lattice Boltzmann method. For both, we explain the algorithm for a bulk fluid and discuss necessary steps to correctly describe inertial focusing.
3.1. Computational fluid dynamics
The numerical treatment of the Navier-Stokes equations has a long history, since only few analytic results are available. The nonlinearities of the underlying equations and the complex, possibly moving, boundaries pose a great challenge. Over the years, methods with different areas of application and with different properties have been proposed. Es-pecially for flows at low Reynolds number, a number of highly specialized and, therefore, efficient methods are available. Examples are the multipole expansion [99], the boundary element method [100], and the use of mobilities in Stokesian dynamics [101]. However, all these methods rely on the linearity of the Stokes equation and are no longer appli-cable for inertial microfluidics. Hence, we focus on general Navier-Stokes solvers which account for the full nonlinearity. In the following discussion, we do not order methods chronologically, but rather by similarity.
3.1.1. Direct discretization of the Navier-Stokes equations
A huge class of numerical methods start from the Navier-Stokes equations (2.8) and discretize the velocity field and the differential operators.
The finite-difference approximation is a particularly popular method to discretize par-tial differenpar-tial equations, such as the Navier-Stokes equations [102]. It describes the solution on a regular lattice and discretizes the differential operators on this lattice using differences of function values between different lattice positions. While conceptually sim-ple, ensuring stability and convergence can be difficult. Furthermore, the regular lattice structure makes it difficult to include boundaries which do not conform to this lattice.