**4. Dynamic clustering and chemotactic collapse in a system of active colloids 49**

**4.7. Relation to the Keller-Segel model**

C(jτ) = 1 n

n

i=1

[q_{6}(t_{i}+jτ)− ⟨q_{6}⟩][q_{6}(t_{i})− ⟨q_{6}⟩]

(∆q_{6})^{2} . (4.12)

Here, {t1, . . . , tn} is a set of equally spaced time points from the stationary state withn
typically around 10000, τ =t_{i}−ti−1, andj ranges from 1 to 1000. We perform a discrete
Fourier transform,

Q_{6}(ω) =

k

j=1

C(jτ) exp(−iωjτ), (4.13)

which, according to Wiener-Khinchine’s theorem, is equal to the power spectrum of q_{6}.
The results for different ζ_{tr} are plotted in Fig. 4.11(a). We fit the spectrum with a
non-normalized Gaussian function and detect its maximum at the position ω_{max}. In the
fluctuating-cluster regime no peak in the power spectrum can be distinguished, rather
more, the curve for ζ_{tr} =−3.2 decreases monotonically. By contrast, in the
oscillating-cluster regime a clear maximum at non-zero frequency ω_{max} exists. We identify the
oscillation state in the state diagram if ω_{max} > 0.01. This value is slightly larger than
zero, in order to being able to clearly identify a maximum. In Ref. [40] the authors
formulated continuum equations for diffusiophoretically coupled active colloids. In the
case where the diffusiophoretic translational velocity acts repulsively, i.e., for ζ_{tr} < 0,
they predict an instability with the onset of spontaneous oscillations. We have shown
that oscillations persist in steady state in a certain region in the state diagram.

As it turns out, activity, rather than phoresis, determines the frequency of the pulsating
cluster. In Fig. 4.11(b) we plot ω_{max} versus ζ_{tr}. From above (fluctuating cluster) and
from below (collapsed cloud) a sharp increase ofω_{max}indicates the onset of the oscillation
regime. The curves for ω_{max} display a plateau-like maximum with a value essentially
independent of ζ_{rot}. For example, the curves in Fig. 4.11(b) belong to Pe = 19 and we
find ω_{max}≈0.012±0.002 for the maximum value. Indeed, the oscillation frequencies are
strongly determined by the activity of the particles. In the inset of Fig. 4.11(b) we plot
the maximum value of ω_{max} versus Pe. Beyond the regime where thermal fluctuations
dominate, which is set by a defined threshold value Pe≈20, the characteristic frequency
exhibits a nearly linear increase in Pe. So, the oscillations become faster if the active
colloids are faster.

We note, that upon increasing Pe the phoretic strengths need to be adjusted as well to
reach the oscillation state. This common scaling of the parameters Pe, ζ_{rot} and ζ_{tr} was
rationalized in Sec. 4.5.3.

Figure 4.11.: (a) The power spectrumQ_{6}(ω) of the bond orientational parameter plotted
against frequency ω (symbols) for different ζtr at ζrot = 4.5. The lines are fits with a
non-normalized Gaussian distribution in order to determine the characteristic frequency
ω_{max}. (b) Characteristic frequency ω_{max} of the oscillating cluster plotted against ζ_{tr}.
Inset: ωmax versus P´eclet number Pe in a double logarithmic plot. The error bars show
the standard deviation of ω_{max} of the maximum values for each Pe.

that such a system can be mapped to the Keller-Segel setup assuming a typical mean field assumption, according to which the two-body position probabilities factorize. However, in our setup we also consider activity and rotational phoresis, which leads to an additional SDE 4.4 for the direction angle φ. In this section, we show how our Langevin equations transform to the celebrated Keller-Segel model by means of a moment expansion. Then, we discuss the implications for our colloidal system.

### 4.7.1. Mapping of the colloidal model to the Keller-Segel equations

First, we transform the equations 4.3 and 4.5 to their assigned Fokker-Planck equation:

∂_{t}P(e,r, t) =−v_{0}∇·(Pe) +ζ_{tr}∇·(P∇c) +D_{tr}∇^{2}P

+ζrot∂φ[(∂φe)·∇cP] +Drot∂_{φ}^{2}P (4.14)
We then derive dynamic equations for the colloidal density P_{0}(r,t) =

P(e,r, t)dφ, as
well as the polar [P_{1}(r, t) =

eP dφ] and the nematic [P_{2}(r, t) =

(e⊗e−^{1}_{2})P dφ] order
parameter.

Integrating Eq. (4.14) over φleads to

∂_{t}P_{0}(r, t) = −v_{0}∇·P_{1}+ζ_{tr}∇·(P_{0}∇c) +D_{tr}∇^{2}P_{0} (4.15)
The term with the polarization P_{1} will ultimately lead to renormalized chemotactic and
diffusion coefficients. To derive a dynamic equation for P_{1}, we evaluate

eEq.(4.14)dφ

and obtain

∂_{t}P_{1}(r, t) =−v_{0}∇P_{2}− v_{0}

2∇P_{0}+ζ_{tr}∇(P_{1}⊗∇c)
+D_{tr}∇^{2}P_{1}−D_{rot}P_{1}+ζ_{rot}(P_{2}− 1

2P_{0}1)∇c (4.16)
where we used ∂_{φ}^{2}e = −e to simplify the rotational diffusion term and ∂_{φ}e ⊗∂_{φ}e =
1−e⊗efor the term withζ_{rot} as coefficient. As in the derivation of Sec. 2.2.3 we neglect
gradients of P1 and P2. However, in our case rotational diffusiophoresis contributes a
term proportional to P_{2}. There is no such term in the expansions of Ref. [37] and [9],
which is why the dynamic equation of the nematic order parameter can be neglected
there. In Ref. [112] it is considered, and in fact our expansion is very much inspired by
Ref. [112]. So, we determine the dynamic equation for the nematic order parameter using

(e⊗e−^{1}_{2}) Eq.(4.14)dφ. Concentrating on the relevant terms, and using the definition
P3(r, t) =

(e⊗e⊗e− ^{3}_{4}e1)P dφ. where e1:= ^{1}_{3}(eαδβγ +eβδαγ+eγδαβ), we obtain:

∂_{t}P_{2}(r, t) =∇(...) +ζ_{tr}P_{2}∇^{2}c

−4DrotP2+ζrot(2P3−3

2P11+1⊗P1)∇c. (4.17)
Now, we close the hierarchy of moment equations by setting P_{3} = 0. Furthermore, as in
Sec. 2.2.3 we neglect the time derivative of the first and second moments on time scales
much larger than the rotational diffusion time 1/D_{rot}, discard gradients of P_{1} and of P_{2}
in sufficiently dilute suspensions, and also neglect higher-order gradient terms inc. Then
we obtain from Eq. (4.17):

P_{2} = ζ_{rot}

8D_{rot}(−3P_{1}1+ 21⊗P_{1})∇c.

We insert this relation into Eq. (4.16). Again, the time derivative is assumed to vanish,
gradients of P_{1} and P_{2} are discarded, and we obtain:

P_{1}(r, t) = 1
2D_{rot}

−v_{0}∇P_{0}−ζ_{rot}P_{0}∇c

+

2ζtr∇^{2}c− ζ_{rot}

4D_{rot}|∇c|^{2}
P1

(4.18)
We concentrate here on perturbations around the uniform and isotropic state and
there-fore neglect in Eq. (4.18) all second-order terms inP_{1} and∇c. Then the polarizationP_{1}
just depends on P_{0}. We substitute it in Eq. (4.15) to obtain the final equation familiar
from the Keller-Segel model:

∂tP0 =ζeff∇(P0∇c) +Deff∇^{2}P0. (4.19)
Here, ζ_{eff} =ζ_{tr}+^{ζ}_{2D}^{rot}^{v}^{0}

rot and D_{eff} =D_{tr}+ _{2D}^{v}^{0}^{2}

rot.

is bounded, γ is irrelevant for the collapse behavior.

In two dimensions, the Keller-Segel model exhibits an instability of the uniform state
towards a chemotactic collapse when its parameters satisfy _{D}^{ζ}^{eff}^{kσ}

cD_{eff} > b, where b is a
con-stant which depends on the geometry of the system (Sec. 3.1). In our unitless parameters
the condition becomes

8πσ(ζ_{tr}+ζ_{rot}Pe)

1 + 2Pe^{2} > b . (4.20)

Judging from this relation, the transition line between gaslike and the collapsed state for
constant Pe parametrized by the parameters ζ_{rot} and ζ_{tr} should be linear. However, for
very large |ζ_{tr}|, the translational phoretic attraction or repulsion between two particles
is stronger than the active motion. In such a situation, either particles collapse and
never detach again (ζ_{tr} ≫0) or they never come close enough to build clusters (ζ_{tr}≪0)
irrespective of ζ_{rot}. At the top left and bottom right of the state diagram Fig. 4.2
one clearly distinguishes the initiation of saturation of the transition line. Hence, the
prediction of linear dependency can only be valid for not too large|ζ_{tr}|and this is indeed
approximately the case in Fig. 4.2.

Furthermore, the condition 4.20 does not reproduce the straight transition line in Fig.

4.8. The reason is that atζ_{rot} <0 pronounced clustering occurs before the collapse takes
place. When testing the transition for ζ_{rot} > 0, where dynamic clustering is absent, we
indeed find the predicted scaling behavior ζ_{tr} ∼ Pe^{2} and ζ_{rot} ∼ Pe of the transition line
reproduced. For this purpose, we performed simulations with positiveζ_{rot}and determined
the parameters where the collapse occured. In Fig. 4.12 we plot ζ_{tr} versus Pe^{2} forζ_{rot} =
1.13. The straight line confirms the expected quadratic behavior. For the inset, we set
ζ_{tr}= 0 and indeed reproduce the predicted linear dependence between ζ_{rot} and Pe.

After all, dynamic clustering cannot be described by the Keller-Segel equation.
Clus-tering occurs when a colloid moves towards the rim of a cluster, then rotates away from it
(ζ_{rot} <0), and finally a delicate balance between activity and diffusiophoretic attraction
sets in: Pe ∼ ζ_{tr}|∇c|. Note that this mechanism crucially depends on ballistic
swim-ming with velocityv_{0}. It cannot be described by our Keller-Segel equation since v_{0} only
appears in effective chemotactic and diffusion parameters, which occur based
approxi-mations assuming very long time-scales. Particularly, desintegration of an agglomeration

æ æ

æ

æ

æ

æ

0 50 100 150 200 250 300 350

0 2 4 6 8

Pe

0 5 10 15 20

0.0 0.5 1.0 1.5 2.0 2.5 3.0

*ζ*rot

Pe^{2}
*ζ*_{tr}

Figure 4.12.: Parameters for which the chemotactic collapse occurs as determined from
simulations. Main plot: ζ_{tr} versus Pe^{2} atζ_{rot} = 1.13. Inset: ζ_{rot} versus Pe at ζ_{tr}= 0
occurs only by diffusion, which is controlled by the diffusion constant D_{eff} ∼ Pe^{2}. This
is the reason that the theoretically predicted collapse line is quadratic in Pe. Assuming
that the colloid attaches directly to two particles at the rim of a cluster, we estimate in
our reduced units |∇c| ≈ 1/(2a/2.33a)^{2} = 1/0.73 and obtain a condition for dynamic
clustering: ζ_{tr} ≈ 0.73Pe. This simple estimate reproduces the straight transition line
between clustering states 1 and 2 in Fig. 4.8 and the region where clustering occurs.