Guo-Jun
Liao
*a,
Carol K.
Hall
b and
Sabine H. L.
Klapp
*a
aInstitut für Theoretische Physik, Technische Universität Berlin, Hardenbergstr. 36, D-10623 Berlin, Germany. E-mail: guo-jun.liao@campus.tu-berlin.de; klapp@physik.tu-berlin.de
bDepartment of Chemical & Biomolecular Engineering, North Carolina State University, Raleigh, NC 27695, USA
First published on 24th February 2020
Based on Brownian Dynamics (BD) simulations, we study the dynamical self-assembly of active Brownian particles with dipole–dipole interactions, stemming from a permanent point dipole at the particle center. The propulsion direction of each particle is chosen to be parallel to its dipole moment. We explore a wide range of motilities and dipolar coupling strengths and characterize the corresponding behavior based on several order parameters. At low densities and low motilities, the most important structural phenomenon is the aggregation of the dipolar particles into chains. Upon increasing the particle motility, these chain-like structures break, and the system transforms into a weakly correlated isotropic fluid. At high densities, we observe that the motility-induced phase separation is strongly suppressed by the dipolar coupling. Once the dipolar coupling dominates the thermal energy, the phase separation disappears, and the system rather displays a flocking state, where particles form giant clusters and move collective along one direction. We provide arguments for the emergence of the flocking behavior, which is absent in the passive dipolar system.
In most studies of active systems with anisotropic interactions, they are assumed to be of short range.11 Less effort has been spent on active systems involving long-range anisotropic interactions, such as dipole–dipole interactions stemming from intrinsic (permanent) dipole moments or dipoles induced by an electric or magnetic field.12 A notable feature of dipolar interactions is that they depend not only on the orientations of the two involved particles, but also on their configuration in space. Specifically, two particles with point dipoles at their center tend to align head-to-tail, where the arrowhead of one dipole moment is directed toward the arrowtail of the other dipole. In contrast, particles configured side by side favor antiparallel alignment. Due to this complexity, one would expect a system of dipolar active particles to show macroscopic structures significantly different from active models with simple polar interactions, such as the Vicsek model8 and its variants.
There are some recent studies, which have focused particularly on the dynamics of dipolar active colloids. For example, ref. 13 and 14 have studied the structural transformation of small dipolar clusters under the impact of activity starting from different (meta)stable configurations, either neglecting13 or including14 hydrodynamic interactions. Another example is a system of spherical, ferromagnetic rollers confined to a fluctuating surface. Here it has been shown, both, experimentally and in simulations, that these particles can exhibit swarming or vortex patterns when energized by a vertical alternating field.15 Owing to the interplay between self-propulsion and dipolar interactions, these vortices persist even on a flat surface.16 In nature, magnetotactic bacteria are known to sense the earth's magnetic field and move along or against the field direction.17–19 Experiments have observed that these bacteria, when confined in a microfluidic channel and placed under an external magnetic field, display clustering behavior along the channel.20 The underlying mechanism of this clustering instability was also investigated analytically and in simulations.21 We also mention dipole-like, active Janus particles with two screened electric charges in each hemisphere. It is found that these particles can self-organize into fingerprint-like patterns at high densities.22 Moreover, an external electric or magnetic field can induce two point dipoles in the respective hemispheres of a Janus particle and thereby further complicate the interactions between particles. By fine-tuning the interactions, it has been reported that the external field can be used to control the collective behavior of Janus particles.23,24 However, these latter models involve even more complex anisotropic interactions as compared to the purely dipolar case. Indeed, the collective behavior of large ensembles of the simplest dipolar active model, that is, dipolar active Brownian particles, is not yet explored. This is the motivation for the present work.
Specifically, we present Brownian Dynamics simulation results for the dynamical self-assembly of a two-dimensional system of active particles, where each particle has a permanent point dipole moment oriented in the plane. In addition, each particle is subject to a self-propulsion force, which is directed along the dipole, as well as to thermal noise. We investigate the collective behavior of our model system for three density regimes, considering a wide range of motilities and dipolar coupling strengths. Our model can be realized, e.g., by Janus particles with a magnetic material coated on one of the hemispheres.25,26 In that case, both the dipole moment and the propulsion force of the Janus particle are directed along its symmetric axis.
For all three densities we present state diagrams illustrating the complex interplay between essentially three phenomena: chain formation (which already occurs in the passive case), motility-induced phase separation, and polar ordering.
The remainder of this paper is organized as follows. In Section 2 we present our model of dipolar active particles and the simulation details, as well as the target quantities and order parameters investigated. Based on analysis of these quantities for a range of parameters, we discuss the collective behavior of the model system in Section 3. Finally, in Section 4 we summarize our findings.
The pair potential between particles located at positions ri and rj (i ≠ j) is of the form
upair(rij,μi,μj) = usr(rij) + udd(rij,μi,μj), | (1) |
![]() | (2) |
In this study, we set the length unit to be σ and fix the repulsive strength ε* = βε = 10, where the thermal energy, β−1 = kBT, is set to be the energy unit (with kB being Boltzmann's constant and T being the temperature). This potential is truncated at the cut-off radius rc = 21/6σ, such that eqn (2) and its derivative vanish to zero continuously at the truncation point. The last term in eqn (1) represents the (long-range) dipole–dipole interaction. Its functional form is given by
![]() | (3) |
We define the dipolar coupling strength as λ = βμ2σ−3 with μ being the magnitude of a dipole moment, i.e., μi = μi.
ṙi = β![]() | (4) |
![]() | (5) |
![]() | (6) |
In eqn (4), the effective propulsion force, which drives the active motion, is given by Fi = F0êi with êi = i. For simplicity, in the remainder of this work we present the impact of the effective propulsion force via the motility v0 = βDtF0. Finally, the random force ξi(t) and torque Γi(t) for the ith particle are zero-mean Gaussian white noise, which satisfy
and
. The angle brackets 〈⋯〉 denote ensemble average, and the symbol ⊗ represents dyadic product.
Eqn (4) and (5) are solved via the Euler–Maruyama method37 with the discrete time step Δt = 2 × 10−5τ, where τ is the Brownian diffusion time, given by τ = σ2/Dt. We choose a quadratic box with side length Lx = Ly = L, and we use periodic boundary conditions in both directions. In order to treat the long-range dipole–dipole interactions, we employ a two-dimensional (2D) Ewald summation, as outlined in Appendix A.1. The particle number is set to N = 1156, unless otherwise stated. All simulations are started with randomly oriented particles located on a square lattice. A typical run then consists of at least 5 × 105 steps for reaching a steady state, followed by a production period of 5 × 105 steps. The statistical properties of the dipolar particles (see Section 2.3) are measured every 500 steps. The simulations are carried out at three values of the mean area fraction Φ = Nπσeff2/(4L2), where the area of a single particle is defined as πσeff2/4. Specifically, we consider the values Φ = 0.12, Φ = 0.23, and Φ = 0.58 (for the exact values, see note ref. 38). To investigate the impact of activity, we perform simulations at different dimensionless propulsion speeds v0* = v0σ/Dt and various dipolar coupling strengths λ. We note that v0* is indeed the same as the commonly used Péclet number Pe.3,39,40
ϕc = 〈nlcl〉/N, | (7) |
The behavior of the global orientational ordering has also been investigated in models of self-propelled particles with “velocity-alignment” interactions, such as the Vicsek model8 and its variants (see, e.g., ref. 9). In these systems, one observes so-called flocking states, where particles gather together and move collectively toward a certain direction. In the present system, the orientational (dipolar) interaction is more complicated (than the Heisenberg-like interactions in the Vicsek model); still, one could imagine that activity-induced flocking states also emerge in the dipolar active system. To this end we consider (as is common in the literature on flocking states9) the parameter
![]() | (8) |
To characterize chain formation in the dipolar active system, we use different strategies depending on the density regime considered. The first strategy is adequate for low densities (Φ ≲ 0.23), where the chains can essentially be considered as isolated objects. In this situation, we consider a chain as a set of, at least, three particles which are mutually “bonded.” Specifically, two particles i, j are regarded as “bonded” in the chain, if the following criteria are fulfilled: |rij| ≤ rp, i·
j > 0, and (
i·rij)(
j·rij) > 0. Here, rp = 1.25σ is set to a distance between the location of the first peak and the first valley of the pair correlation function.58–60 Based on these rules, we quantify the low-density chain formation via the parameter
ϕp = 〈Np〉/N, | (9) |
At high densities (e.g., Φ = 0.58), the particles are obviously closer to each other, and the string-like structures are no longer isolated. This causes the degree of polymerization ϕp to no longer be appropriate to characterize the self-assembly. To overcome this difficulty, we propose a different order parameter to describe the chain formation quantitatively. The starting point is the dipole–dipole correlation function
![]() | (10) |
![]() | (11) |
Here, the distance rs is set to rs = L/8 ≈ 5σ such that rs is much smaller than the half of the box size (which corresponds to the length scale of the periodic boundary conditions). In the presence of chain-like structures with head-to-tail alignment of neighboring particles, we expect the angular correlation function, μ(θ), to display positive maxima both, in front of and behind the reference particle (i.e., at θ = 0 or π). In contrast, the dipole moments of the particles on the right- and left-hand side of the reference particle should remain rather uncorrelated (i.e.,
μ(θ) ≈ 0 for θ = π/2 or −π/2). With this picture in mind, we measure the quantity
Z = ![]() ![]() | (12) |
State | Clustering | Orientational ordering | Chain formation | |
---|---|---|---|---|
Φ ≲ 0.23 | Φ = 0.58 | |||
Homogeneous, isotropic fluids | ϕ c ≤ 0.5 | ϕ e ≤ 0.5 | ϕ p ≤ 0.5 | Z ≤ 0.17 |
Chain-like structures | ϕ c ≤ 0.5 | ϕ e ≤ 0.5 | ϕ p > 0.5 | Z > 0.17 |
Micro-flocking | ϕ c ≤ 0.5 | ϕ e > 0.5 | — | |
Macro-flocking | ϕ c > 0.5 | ϕ e > 0.5 | — | |
Motility-induced clustering | ϕ c > 0.5 | ϕ e ≤ 0.5 | — |
![]() | ||
Fig. 2 Representative simulation snapshots at Φ = 0.12. Particles are colored according to their dipole orientations as indicated by the color ring in the inset. |
To quantitatively characterize the motility-induced destruction of chain-like structures, we plot in Fig. 3 the degree of polymerization ϕp as a function of the motility v0* for various coupling strengths λ. In the absence of dipole–dipole interactions (λ = 0), particles do not self-assembly into chain-like structures, still we observe a slight increase in ϕp [see eqn (9)] upon increasing motilities v0*. This can be attributed to the dynamical clustering of active Brownian particles into small (finite-sized) aggregates.3,36,61 Indeed, if the particles inside these finite-sized clusters fulfill our criteria of being bonded, the order parameter ϕp will be non-zero (yet small), despite the fact that dipole–dipole interactions are absent. Once the dipole–dipole interactions are introduced, pronounced chain-like structures appear for coupling strengths λ ≳ 6 and low motilities v0*, giving rise to large values of ϕp (ϕp ≈ 1). Upon an increase in the motility, the degree of polymerization ϕp gradually decreases, reflecting that the self-propulsion opposes the formation of chain-like structures. This behavior resembles that seen in passive, dilute systems of self-assembled dipolar particles, when the temperature is increased.62 In this sense, the particle motility in the active dipolar system may be viewed as an analog to the temperature in the passive system.
![]() | ||
Fig. 3 Degree of polymerization ϕp as a function of the motility v0* at Φ = 0.12 for the coupling strength λ = 0 (black circles), 4 (red squares), 6 (green diamonds), and 10 (blue triangles). The dashed horizontal line at ϕp = 0.5 marks our criterion for a state with chain-like structures (see Table 1). The solid lines are guides to the eye. |
The degree of polymerization ϕp only describes the chaining behavior globally. To supplement our analysis of chain formation, we measure the chain size distribution function P(n), which represents the probability that a randomly selected chain consists of n particles. In the presence of chain formation, the simulated systems usually contain several long chains and many short chains, such that P(n) becomes very small at large n. Therefore, we plot in Fig. 4 the chain size distribution weighted by n, where nP(n) is proportional to the probability that a randomly selected particle belongs to a chain with a size n. As seen in Fig. 4(a), the weighted distribution curves of non-dipolar particles vanish at around a chain size n ≲ 10, and the values slightly increase with increasing motilities. This is consistent with the aforementioned dynamical clustering. The corresponding functions at large dipolar coupling strengths (λ = 10) look very different, as shown in Fig. 4(b). Here, we observe a broad peak between 10 ≲ n ≲ 100 at v0* = 0, reflecting formation of long chains. Upon an increase in the motility, this peak is gradually shifted to a smaller chain size. Finally, this peak disappears once v0* ≳ 60. This suggests a vanishing of chains with the most probable size at large values of v0*.
Inspecting again Fig. 2(d) we see that even at the largest coupling strength considered, the dipole moments do not align on a length scale comparable to that of the simulation box, regardless of the values of v0*. In other words, there is no pronounced global orientational order. Indeed, through measuring the orientational order parameter ϕe, we confirm that ϕe remains small (ϕe < 0.5) for the parameter combinations explored in Fig. 1, indicating that the systems are globally isotropic at low densities.
![]() | ||
Fig. 6 Representative simulation snapshots at Φ = 0.23. Particles are colored according to their orientations as indicated by the color ring in the inset. |
We now turn to the impact of the motility. At small values of λ, an increase in v0* toward high values leads to the formation of finite-sized clusters [see Fig. 6(b)] which are, however, not stable. Thus, the system remains homogeneous and isotropic on average. This finding is consistent with earlier research on active Brownian particles. In particular, the phenomenon of motility-induced phase separation, accompanied by formation of “giant” clusters (ϕc ≈ 1), is known to occur only at higher densities.3 As Fig. 5 reveals, the picture changes if the dipolar coupling becomes larger. Consider, e.g., the case λ = 6. Increasing the motility from zero, the chain-like structures observed at small v0* break [see Fig. 6(c)], and the particles instead form small clusters characterized by the same orientation [see Fig. 6(d)]. More quantitatively, the clustering parameter ϕc remains small (ϕc < 0.5) but the orientational order parameter ϕe reaches values above 0.5.
The full behavior of the orientational order parameter ϕe as function of v0* is shown in Fig. 7. For strongly coupled systems (λ ≥ 3), the order parameter abruptly changes from (essentially) zero to large values at a “critical” motility v0,c*(λ) ≈ 20. We take this as an indication for a motility-induced formation of ordered, yet small, clusters. According to Table 1, this is characteristic of a “micro-flocking” state. We recall that, for active particles, the development of average alignment between neighboring particles implies that they move along the same direction. To explore whether these small clusters are indeed “microscopic” structures, we measured the cluster size distribution. It turns out that the characteristic cluster size n0 does not scale up with the system size N (see the detailed analysis in Appendix A.2).
![]() | ||
Fig. 7 Global polarization ϕe as a function of the motility v0* at Φ = 0.23 for the coupling strength λ = 0 (black circles), 2 (red squares), 3 (green diamonds), 4 (blue triangles up), and 6 (orange triangles down). The dashed line indicates the value corresponding to the emergence of flocking states (see Table 1). The solid lines are guides to the eye. |
Inspecting Fig. 7 again, it is found that for systems at an intermediate coupling strength λ = 2, the order parameter ϕe gradually increases from zero as v0* ≳ 100, suggesting that the “micro-flocking” state might still appear. However, further increase in v0* requires to employ even smaller simulation time steps Δt < 2 × 10−5τ, so as to prevent the numerical instability. Combining this aspect with the time-consuming Ewald summation, it turns out that the simulations become computationally unfeasible for investigating the flocking state transition at small dipolar coupling strengths (λ < 3).
A behavior similar to the “micro-flocking” state has been previously observed in systems of active Brownian disks with polar alignment.9 Indeed, the density Φ = 0.23 considered here is close to the effective packing fraction ϕ* = 0.256 defined in ref. 9, where a comparable state with ‘microscopic’ polar clusters has been discovered at high motilities and intermediate polar coupling. Upon further increase in polar coupling, the system in ref. 9 displays macroscopic structures, such as moving patterns of bands or lanes. Here, we did not find such patterns in the parameter range explored. In this context, we also mention a state observed in the Vicsek model at low densities and low noise (yielding strong alignment).8
In this state, the point-like particles form small groups within which the particles move together along one (random) direction. However, the system does not show global ordering. In both of these models, i.e., active polar disks and the Vicsek model, the particles favor alignment of the propulsion direction with their neighbors, independent of their relative positions. This is different from the particles considered here whose interaction depends not only on the relative alignment but also on their configuration.
In order to gain a deeper insight into the flocking behavior of the present system, we draw in Fig. 8(a–d) sketches of the motion of two highly motile, strongly coupled particles starting from four representative configurations. The initial configurations (related to time t0) in Fig. 8(a) and (b) are those with the lowest dipolar interaction energy. These configurations are therefore stable (up to thermal fluctuations) in the passive case. Upon switching on the motility (v0* > 0), the particles in the head-to-tail configuration [see Fig. 8(a)] will move along a straight line in the same direction (as already argued in ref. 13), until thermal fluctuations set in. We thus consider this configuration as relatively stable not only in the passive, but also in the active case. This is different in the antiparallel side-by-side configuration depicted in Fig. 8(b) at t0. Here, activity leads to motion in opposite direction, yielding this configuration unstable. In Fig. 8(c), the initial configuration at t0 is energetically unstable. Still, at least for short time, one expects the particles to move together along the same direction until the dipolar repulsion drives them away from one another. Finally, in Fig. 8(d) the initial configuration is a slightly distorted (induced by rotational noise) head-to-head configuration. This leads to a strong dipolar torque which, combined with strong activity, pushes the particles apart. Combining these arguments, we conclude that the most stable configuration for two dipolar active particles is the head-to-tail alignment. While this is similar to the passive case, a major difference is that the antiparallel side-by-side configuration [Fig. 8(b)] breaks immediately apart when the particles become active, yielding this configuration unstable. At low densities, we therefore expect to find short straight chains moving along their individual long axis, consistent with the observation in the Fig. 2(d). Upon an increase in the density, short chains collide with each other more frequently than at low densities. As explained above, antiparallel alignment of the short chains is not stable. Therefore, once two short chains collide, they tend to align along the same direction and, thus, form a dynamical polar cluster, as illustrated in Fig. 8(e). At a sufficiently high density (e.g., Φ ≳ 0.23), the frequent collisions between short chains finally lead to a “micro-flocking” state, in which the dipolar active particles form finite-sized clusters with a polar order [see Fig. 6(d)]. This is fundamentally different from the alignment mechanism in other models such as the ones with ferromagnetic-like interactions, where the configuration at t0 in Fig. 8(c) is kinetically stable.9,10
![]() | ||
Fig. 10 Representative simulation snapshots at Φ = 0.58. Particles are colored according to their orientations as indicated by the color ring in the inset. |
The motility-induced phase separation can be explained by a self-trapping mechanism: When a highly motile particle enters into a dense region occupied with other particles, this particle is temporally slowed down due to frequent collisions. Such a slowing-down effect makes this region even denser, thus creating a positive feedback loop which leads to the formation of giant clusters.39 We note that the particle orientations in the two phases are essentially uncorrelated for non-dipolar particles (λ = 0), as seen in Fig. 10(b).
To get a first impression of the impact of dipole–dipole interactions on the clustering behavior, we plot in Fig. 11 the fraction of the largest cluster, ϕc [see eqn (7)], as a function of the motility v0* for various coupling strengths λ. Consistent with the motility-induced phase separation, the curves for non-dipolar and weakly coupled active particles (λ = 0–1) display a sharp increase from zero at v0* ≈ 20 and reach large values (ϕc ≈ 0.8) at v0* = 100. The curves reach ϕc = 0.5 at larger v0* with increasing λ, which is also reflected by the boundary between black dots and red squares in Fig. 9. For strong dipolar coupling (λ ≳ 2), the values of ϕc at high motilities remain in the same range (ϕc ≈ 0.8, indicating again formation of large clusters), but the increase in ϕc at small v0* is much less pronounced. This already indicates a strong impact of dipolar interactions. For more detailed analysis of the clustering behavior, we discuss the cluster size distribution in Appendix A.2.
To gain further information, we calculate position-resolved local area fractions, ϕ, based on a Voronoi tessellation.36,40 (We note that, contrary to ref. 36, we did not perform a short time average of ϕ, since the polar clusters characterizing the dense state migrate over time and, therefore, are not stationary within the short-time interval. The flocking behavior will be discussed later in detail in Section 3.3.2.) Fig. 12 shows the probability distribution of the local area fractions, P(ϕ), for different coupling strengths at v0* = 100. For non-dipolar active particles (λ = 0), P(ϕ) reveals a clear double-peak structure that reflects the separation between the dilute and the dense phase. The coexisting densities correspond to the location of the two local maxima of P(ϕ). Upon an increase in the coupling strength up to λ = 2, the double-peak structure gradually disappears. Instead, we observe the emergence of a single peak located at a density slightly larger than the mean density (Φ = 0.58), as well as a broad shoulder on the left. This suggests the disappearance of the phase separation observed at zero and small λ.
The same conclusion can be drawn from Fig. 13, where we plot the coexisting densities in the (v0*, ϕ) plane for various λ. The location of the high-density branch at λ = 0 can be explained as follows: At high motilities v0*, the non-dipolar system (λ = 0) exhibits “giant” clusters, composed of randomly oriented particles [see Fig. 10(b)]. These particles are separated by a distance close to the effective hard sphere diameter σeff defined in Section 2.2. As a result, the local area fraction in the densely packed region should be close to the close-packing fraction, . Upon a slight increase in dipolar coupling, the coexistence branches for the dilute region, ϕgas, are shifted toward higher densities, while the branches for the dense region, ϕden, move toward lower area fractions. As a result, the area surrounded by the curves of ϕden and ϕgas in Fig. 13 significantly shrinks with increasing λ, indicating that motility-induced phase separation is generally suppressed by the dipolar interactions. Within the present simulations, the phase separation and the corresponding coexistence curves disappear once λ > 1. We note that this observation is in contrast with the findings in a recent study of active Brownian particles with additional interactions of velocity alignment at high densities.10 This shows that different types of orientational interactions may have entirely different impacts on motility-induced phase separation.
Yet another perspective on the disappearance of the phase separation due to dipolar interactions emerges when we consider the self-trapping mechanism (which plays a key role for clustering at λ = 0). To this end, we compute the normalized speed of the ith particle, vi*/v0*, versus the local density ϕi = ϕ(ri), where the speed is given by vi* = |(Δri/σ)/(Δts/τ)| with Δts = 10−2τ. After averaging over N particles and over at least 1000 snapshots, we plot in Fig. 14 the normalized particle speed versus the local density. For non-coupled (λ = 0) and weakly coupled (0 < λ ≲ 1) systems, we observe a linear decay,
v*(ϕ)/v0* = 1 − adϕ, | (13) |
![]() | ||
Fig. 14 Normalized particle speed v*(ϕ)/v0* as a function of the local area fraction ϕ for the coupling strength λ = 0 (black circles), 0.5 (red squares), 1 (green diamonds), 2 (blue triangles up), and 5 (orange triangles down) with the motility v0* = 100. The solid lines represent fits to eqn (13). Inset: Fitted decay amplitude ad as a function of λ. The line is drawn as a guide to the eye. |
To obtain a more complete (yet qualitative) picture on the role of the function ad(λ) for the appearance of the motility-induced phase separation, we consider the effective free energy proposed in ref. 41,
f(ϕ) = f0(ϕ) + frep(ϕ), | (14) |
![]() | (15) |
frep(ϕ) = krepΘ(ϕ − ϕt)(ϕ − ϕt)4 | (16) |
For conventional active Brownian particles (λ = 0), the decay amplitude, as a crucial parameter describing the degree of self-trapping, is typically chosen to be ad = 1 for 2D and ad = 1.3 for 3D.64 With this in mind, in Fig. 15 we plot the effective free energy f(ϕ) for the decay amplitude ad = 1.08, 0.95, and 0.8. At low densities, f(ϕ) is nearly independent of ad, while at high densities, f(ϕ) increases with decreasing ad. From f(ϕ), we determine the coexisting densities through the common tangent construction. To compare the influence of the decay amplitude, ad, on the phase separation from the theoretical perspective and the simulation results, we plot in the inset of Fig. 15 the coexisting densities in the (ad,ϕ) plane for the effective free energy (blue circles) and simulations at v0* = 100 and 0 ≤ λ ≤ 1 (orange squares). More specifically, for the above simulations we plot the coexisting densities determined from Fig. 12versus the associated ad plotted in the inset of Fig. 14. The coexisting densities obtained from both methods describe the same trend, that is, a decrease in the decay amplitude ad causes a reduction in ϕden and a growth in ϕgas. Further, the phase separation disappears once v(ϕ) decays sufficiently slow (ad ≲ 0.6 for the effective free energy, and ad ≲ 0.5 for simulations). This clearly indicates that the motility-induced phase separation is suppressed when the function v(ϕ) decreases too slowly, i.e., when the dipolar coupling becomes too strong.
![]() | ||
Fig. 16 Magnitude of the average orientation ϕe as a function of the coupling strength λ at Φ = 0.58 for v0* = 0 (black circles), 20 (red squares), 60 (green diamonds), 80 (blue triangles up), and 100 (orange triangles down). The dashed line ϕe = 0.5 indicates the value above which the particles form a flocking state (according to Table 1). The solid lines are drawn as a guide to the eye. |
To better evaluate the dependence of orientational correlations on the direction of the connecting vector rij, Fig. 17(e) shows the angular correlation function, μ(θ) [defined in eqn (11)], for the strong dipolar coupling strength λ = 5. As mentioned in Section 2.3.3, the difference between the maximum
μ(θmax) and the minimum
μ(θmin) of the curve, Z [defined in eqn (12)], can be used as a measure of chain formation. Upon increasing v0*, Z decreases and, at the same time, the curves in Fig. 17(e) are shifted toward lager positive correlations. An overview of the behavior of Z as a function of v0* (for various λ) is given in Fig. 18. The data curves bear close resemblance to those for ϕp in Fig. 3. For a non-dipolar system (λ = 0), where chains are absent, the order parameter Z does not show any significant behavior within the range of motilities explored. In contrast, at large coupling strengths (λ = 4 and λ = 5), the decrease in Z with increasing v0* reflects the disappearance of a state with chain-like structures. Specifically, the curves for λ ≥ 4 show a sharp decrease, accompanied by a point of inflection Z ≈ 0.17. We choose this as a threshold value. We note that the simulations of nearly passive, strongly-coupled systems (e.g., λ = 5 and 6) are severely plagued by large fluctuations of the order parameter Z over time. In these cases, we were unable to determine whether or not the simulations had reached a steady state even after a very long simulation time (t > 300τ). The corresponding region in the state diagram in Fig. 9 is marked by dashed lines.
![]() | ||
Fig. 17 Spatial correlation function of dipole moments, gμ(r⊥, r‖) defined in eqn (10), at Φ = 0.58 and λ = 5 for the motility v0* = 0 (a), 8 (b), 16 (c), and 100 (d). The black circle represents the reference particle with the red arrow indicating the particle's orientation as well as the dipole moment. The area ![]() ![]() |
When the motility is small, the state diagrams are similar for all densities considered. Specifically, homogeneous and isotropic fluids are observed for nearly passive particles with weak dipolar coupling, whereas strong dipolar coupling leads to chain-like structures. At high motilities, the state diagrams strongly depend on the mean area fraction Φ. At low densities (Φ = 0.12) and strong dipolar coupling, an increase in the motility v0* from zero causes chain-like structures to break into fragments of short chains and individual beads. At intermediate densities, passive, strongly-coupled dipolar particles self-assemble into a state with chain-like structures. With increasing v0* from zero, chains start to break and the system displays a micro-flocking state, where particles form finite-sized clusters with polar order. Finally, at a high densities and strong dipolar interactions, we observe a motility-induced transition from a state with chain-like structures into a micro-flocking state and finally into a macro-flocking state, where particles show global orientational ordering and form “giant” clusters.
To provide a simple argument for the emergence of polar order (which is absent in the passive 2D case), we have considered the time evolution of two dipolar particles starting from four representative initial configurations. As a result of the interplay between dipolar coupling and self-propulsion, the head-to-tail configuration remains the most stable one (same as the passive case), while the antiparallel side-by-side configuration is destabilized. With increasing particle number, the head-to-tail alignment mechanism can finally lead to a flocking state if the density is sufficiently high (Φ ≳ 0.23).
It has been reported earlier that finite size effects may hide crucial features of the flocking behavior in active systems, such as the first order nature of the flocking transition followed by the formation of traveling bands.6 We note that the sizes of the present simulations are limited to the order of 103 particles, due to the long-range character of the dipolar interactions and, subsequently, the expensive computational cost. Besides the qualitative behavior, it would be very interesting to study the scaling behavior and explore whether the dipolar active particles belong to any of the existing universality classes, such as the one of the Vicsek model.4,66,67 To this end, however, it would be necessary to perform extensive simulations with system sizes much larger than N ≈ 103, which is, again, limited by the computational resources. Nevertheless, it is worth mentioning a fundamental difference regarding the flocking transition: in the Vicsek model the onset of the flocking state is characterized by simultaneous spontaneous appearance of density and orientational inhomogeneities.6–8 In contrast, in our model the emergence of global orientational order is decoupled from the formation of large-scale structures, such as “giant” clusters.
At Φ = 0.58 and large values of v0*, the system displays motility-induced phase separation, where “giant” clusters composed of densely-packed particles with random orientations coexist with freely moving particles in the dilute region. The phase separation persists as long as the dipolar coupling λ is negligible against thermal fluctuations (0 < λ ≲ 1). Once the dipolar interactions dominate (λ > 1), the orientations of the active particles are no longer uncorrelated and the particles tend to align with their neighbors, leading to a break down of the self-trapping mechanism and a subsequent suppression of the phase separation.
We note that our model does not account for hydrodynamic interactions between the particles. In the absence of dipolar coupling, simulation studies have reported that hydrodynamic interactions tend to suppress motility-induced phase separation due to, either, the near-field interactions,68,69 or the rapid decorrelation of the particle orientations.70–72 Through investigating the pair distribution function, it has been shown that hydrodynamic interactions generally damp out the translational structure of active particles at high densities.73 Moreover, hydrodynamically interacting particles with a certain range of force dipole strengths can spontaneously form a state with global polar order,68,69,74,75 which can be attributed to either the actively induced rotation-translation coupling,76 or the near-field lubrication forces.68,69 Therefore, one may expect that dipolar coupling and hydrodynamic interactions combine to further suppress the phase separation and promote the polar ordering. The details of the resulting collective behavior are a topic of future studies.
In the real world, active particles are often asymmetric.77,78 As a result, the effective propulsion force does not coincide with the particle's center of mass, thus generating a propulsion torque which induces chiral active motion.35,79 Such a mechanism is also expected for the magnetic or dielectric Janus particles, whose dipole moments (either permanent or induced) are mostly shifted from the center of the individual particles.80–83 Therefore, it would be very interesting to investigate models accounting for the chiral motion and the shifted dipole moment, so as to obtain more comprehensive understanding of dipolar active systems.
![]() | (17) |
![]() | (18) |
![]() | (19) |
![]() | (20) |
The Langevin eqn (4) and (5) involve the forces and torques due to the dipolar interactions. From eqn (17), the force acting on the ith particle due to dipole–dipole interactions is given by
Fi,dd = −∇riUdd = FRi,dd + Fk≠0i,dd, | (21) |
![]() | (22) |
![]() | (23) |
![]() | (24) |
In the present work, every dipole moment is restricted to point in the xy-plane with ψi characterizing its direction relative to the x-axis. Therefore, the torque acting on the ith particle is always directed along the z-axis. It is given by
Ti,dd = −∂ψiUdd![]() | (25) |
![]() | (26) |
![]() | (27) |
nP(n) = P(1)n−1e−(n−1)/n0, | (28) |
![]() | ||
Fig. 19 Weighted distribution of cluster size at Φ = 0.23 for the particle number N = 100 (black circles), 225 (red squares), 400 (green diamonds), 1156 (blue triangles up), and 2500 (orange triangles down). The orange solid line indicates the curve fitted to N = 2500 according to eqn (28). Inset: Characteristic size of clusters, n0, as a function of particle number N. |
To provide additional information on the clustering behavior at high densities (Φ = 0.58), we plot in Fig. 20 the weighted cluster size distribution for non-dipolar and strong coupled dipolar active particles at Φ = 0.58. At the coupling strength λ = 0, our model reduces to the limiting case of active Brownian particles.3,63 Here, the correspondent weighted cluster size distribution function undergoes a transition from an exponential decay (v0* ≲ 20) into a curve with a power law decay at small cluster sizes n and a peak at large n (v0* ≲ 20), as shown in Fig. 20(a). Interestingly, the distribution functions of strongly-coupled dipolar particles (λ = 5) bear resemblance to the case of active Brownian particles, as can be seen in Fig. 20(b). For passive dipolar particles (v0* = 0), the weighted distribution function vanishes at n ≈ 30, which is three times as large as in the non-dipolar case. This is due to the fact that the strongly-coupled particles tend to form head-to-tail configurations and therefore display chain-like structures. With increasing motilities, such structures disappear. Further, the dipolar particles start to align with their neighbors and form into clusters, as shown in Fig. 10(d). As indicated by the peaks at n ≈ 103 in Fig. 20(b), these particles form “giant” clusters when v0* ≳ 20.
This journal is © The Royal Society of Chemistry 2020 |