Xiaoyue
Wu
a,
Katherine
Skipper
b,
Yushi
Yang
b,
Fergus J.
Moore
b,
Fiona C.
Meldrum
a and
C. Patrick
Royall
*c
aSchool of Chemistry, University of Leeds, Woodhouse Lane, Leeds, LS2 9JT, UK
bH. H. Wills Physics Laboratory, University of Bristol, Bristol BS8 1TL, UK
cGulliver UMR CNRS 7083, ESPCI Paris, Université PSL, 75005 Paris, France. E-mail: paddy.royall@espci.fr
First published on 12th March 2025
Colloidal particles self assemble into a wide range of structures under external AC electric fields due to induced dipolar interactions [Yethiraj and Van Blaaderen, Nature, 2003, 421, 513]. As a result of these dipolar interactions, at low volume fraction the system is modulated between a hard-sphere like state (in the case of zero applied field) and a “string fluid” upon application of the field. Using both particle-resolved experiments and computer simulations, we investigate the emergence of the string fluid with a variety of structural measures including two-body and higher-order correlations. We probe the higher-order structure using three-body spatial correlation functions and a many-body approach based on minimum energy clusters of a dipolar-Lennard-Jones system. The latter constitutes a series of geometrically distinct minimum energy clusters upon increasing the strength of the dipolar interaction, which are echoed in the higher-order structure of the colloidal fluids we study here. We find good agreement between experiment and simulation at the two-body level. Higher-order correlations exhibit reasonable agreement between experiment and simulation, again with more discrepancy at higher field strength for three-body correlation functions. At higher field strength, the cluster population in our experiments and simulations is dominated by the minimum energy clusters for all sizes 8 ≤ m ≤ 12.
Colloidal dipolar systems fall broadly into two categories. Some, for example ferromagnetic nanoparticle systems, like atoms and molecules, have an intrinsic dipole moment6 and can be modeled with the Stockmayer model which combines a dipolar interaction with a Lennard-Jones interaction.7 These systems exhibit intriguing string-like structures,8,9 with branching,10 coiling11 and clustering7 behavior, not to mention a ferromagnetic transition.12 Rather than spontaneous dipolar interactions, in other colloidal systems, dipoles may be induced by an external electric or magnetic field.4,5,13 This has the consequence that the dipolar interactions are aligned in the direction of the applied field. Using ferromagnetic and superparamagnetic nanoparticles in an external field then opens further possibilities such as a very strong response to the field.14 Other, more exotic possibilities include the use of a biaxial field, leading to phenomena such as in-plane condensation in (quasi) 2d systems,4,15 control of assembly and tuning of interactions16,17 and direct observation of cluster growth.18 In addition to their fundamental interest, such dipolar colloidal systems may find application as electrorheological fluids,19,20 hydraulic valves21 and photonic materials.22
Here we shall focus on dipolar systems with an external field. A particular attraction of these systems is the ease with which the dipolar interactions can be tuned with the external field. Indeed a combination of real23 and reciprocal space24,25 studies of such systems enables the investigation of a variety of crystal structures including fcc, hcp, bcc, body-centered tetragonal (bct) and body-centered orthorhombic (bco) structures,5,23,26 along with a transient labyrinthine structure.20 Tuning the electric (or magnetic) field in situ enables the control of phenomena such as a martensitic transition.27 Adding softness28 or attractions29 to the interaction potential further increases the range of structures into which the system may self-assemble.
In addition to the rich crystalline phase behavior, dipolar colloids feature a fluid phase at lower colloid volume fraction and electric field strength than those at which the crystals are found. Interestingly, the symmetry-breaking dipolar interactions cause this fluid to assemble into string-like structures which are aligned in the direction of the electric field.20,23,30 This “string fluid” has been investigated analytically31 and can form the basis for producing “colloidal polymers”.32,33 Meanwhile, it is possible to take the system out of equilibrium, which enables investigation of string growth mechanisms,18 and aggregation phenomena between the strings.34 This suggests that the structure of the “string fluid” may be rather interesting and it has been investigated analytically31 and in reciprocal space9 in addition to real space.23,29
In atomic and molecular systems, characterizing structure in the fluid state beyond pair correlations is challenging, although not impossible.35,36 With colloidal systems, particle-resolved studies4 which deliver coordinates in real space are amenable to the measurement of higher-order correlation functions, relevant to a variety of phenomena such as dynamical arrest35,37–39 and polymorph selection40–42 and crystal precursors.41,43 Furthermore, theoretical treatments have been developed to describe the higher-order structure of hard sphere colloids.44–46
Methods to characterize higher-order structure include three-body correlation functions such as g3,44,47,48 and higher-order correlations such as common neighbor analysis (CNA)49 and Voronoi face analysis.50 These have been shown to be successful in studying the structure of fluids and glasses.51,52 Another strategy, the bond orientation order (BOO) parameters developed by Steinhardt et al.,53 focuses on the local symmetry around a central particle. This method has been shown to be very useful in the study of crystallization, especially in identification of small crystalline clusters in a supercooled liquid54,55 and also in the characterization of fivefold symmetric order in amorphous systems.37 With the popularity of machine learning rising, it has been applied effectively to local structure in amorphous materials, for example by combining many structural metrics such as the pair correlation function.41,56 Other examples include combining it with local descriptors such as CNA and BOO to better characterize the local environment around a single particle in disordered materials.41,57,58 Both supervised59 and unsupervised learning60 have been used to further our understanding of supercooled liquid and glass forming systems.
The methods discussed above are geometric in nature. An alternative approach, which takes into account the interactions between the constituent particles of the system, has its roots in the work of Sir Charles Frank,61 who postulated that since the minimum potential energy configuration of 13 Lennard-Jones atoms corresponds to an icosahedron, that this would be a common geometric motif in (supercooled) liquids. With the advent of energy landscape calculations,62 it has become possible to determine the structure of minimum potential energy clusters for a wide range and size of systems including the Lennard-Jones,63 “Sticky Sphere”,64 Stockmayer65 models and more exotic interactions relating to multiple fields.66 Since the dipolar interaction of the Stockmayer model is not constrained to lie in any particular direction it thus corresponds to a molecular (or nanoparticle7) system, rather than a colloidal dipolar system in an external field in the context of the discussion above.
Identifying local arrangements of particles in bulk systems whose bond network is identical to such clusters can be carried out using the topological cluster classification (TCC).67,68 The TCC has been used to identify locally favored structures or minimum energy clusters in systems undergoing dynamic arrest,38 colloid–polymer mixtures interacting via the Morse potential,69,70 colloidal suspensions with attractive interactions,71 colloidal gels48,72 and the liquid–gas interface.73 However, thus far this method has only been applied to systems with a spherically symmetric interaction, and those whose state point is fixed by the composition of the system. Some of us have recently determined the minimum energy clusters for a Lennard-Jones-dipolar system where the dipoles are induced in a particular direction (Fig. 1),68 which opens the possibility to use this method to probe the higher-order structure of dipolar colloids. in this way, we can apply the TCC method to an experimental system with asymmetric interactions which can be tuned in situ.
![]() | ||
Fig. 1 Rigid minimum energy clusters of the dipolar-Lennard-Jones system for various sizes m. 8B, 9B, 10B, 11C and 12B are minimum energy clusters for the Lennard-Jones system (γ = 0). Different geometries correspond to minimum energy clusters as a function of the dipolar strength γ. Here we consider rigid clusters only. The clusters are formed from rings of three, four or five particles. These are coloured grey. In the axis perpendicular to the rings are so-called spindle particles, colored yellow, one above and one below the ring. Single spindle particles are coloured red.67,68 |
Herein, we report a combined experimental and computer simulation study which carries out such a study of higher-order structure of dipolar colloids. Since the electric field can be tuned at will, it is possible to vary the state point of the system in situ. This is somewhat unusual for colloidal systems, where the state point is often fixed by the composition of the system. Here we explore the equilibrium string fluid phase, but we can also increase the field such that the system becomes metastable to a phase coexistence between a fluid and a body-centered tetragonal crystal.26 We consider pair correlations in the form of radial distribution functions g2(r) and three-body correlations in the form of order parameters to determine “string-like” configurations and also the triplet correlation function g3(r,r′,η). We use the topological cluster classification67,68 to explore higher-order spatial correlations in the form of minimum energy clusters of the dipolar-Lennard-Jones interaction.68 Before concluding, we provide an outlook on the study of non-equilibrium behaviour accessed via higher-order correlations.
When colloids are subjected to an external electric field E, a dipole–dipole interaction is induced between the particles which reads
![]() | (1) |
![]() | (2) |
![]() | (3) |
Combining with the hard sphere interaction noted above uhs, the total interaction between two colloids under the electric field becomes
utotal(r, θ) = uhs(r) + udip(r, θ). | (4) |
We determine the dipolar contribution to the interaction potential between the particles by evaluating eqn (2) and (3) with the particle diameter σ, the solvent dielectric constant 83 and the measured value of the local electric field E.
In order to construct the sample cell to hold the colloidal suspension, two indium tin oxide glass slides were separated with spacer silica particles of approximately 60 μm in diameter. This created a transparent sample cell with electrodes top and bottom. The electrodes were connected to a signal generator to provide an AC electric field across the cell.
A Leica SP8 confocal microscope was used to monitor the system under the applied electric field. During each measurement, a stack of 3d confocal images of at least 200 “slices” of xy images along the z direction were taken with 256 × 256 pixels with at least ten pixels per particle diameter in all directions. 3d images were acquired at least 30s apart. Following related work with dipolar colloids,23,83 only particles at least ten diameters from the wall were analyzed, to ensure that there were no significant wall effects. We saw no influence from the wall proximity in any of our measurements and we conclude that we can treat the system as bulk, as is typical for such particle-resolved studies. It is worth noting that the system sizes used in this experimental technique are not huge.4 The shortest dimension of the system (here 60 μm) is two orders of magnitude larger than that of the colloidal particles.
In each measurement, the applied voltage and the thickness of the electrical cell was measured in order to determine the electric field strength. The frequency of the electric field applied was 1 MHz and the field strength was up to 0.3 V μm−1 measured peak-to-peak. Before each measurement, the system is allowed to stabilise for at least 20 min. The Brownian time τB = (3πνσ3)/(4kBT) ≈ 6.09 s for our system, so we consider this time quite sufficient to relax equilibrium states. Here ν is the solvent viscosity. At least 50 3d images each separated by 30s were taken for each state point, in the same place in the sample. In equilibrium, this provides sufficient statistics. However, out of equilibrium, for high field strengths, there will likely be some dependence on the history of the system, to which we return below in Sections 4 and 6.
To mitigate such blurring in the z-direction, we refine the first set of coordinates determined as follows. Knowing the size of the particles, the algorithm predicts an image based on the set of particle coordinates. This predicted image is then iteratively compared with the original measured image and the coordinates moved following a Monte Carlo method using the difference in pixel values between the predicted and measured image to minimize the differences between them. Further details of our method, including the source code may be found in Yang.87
Colloid tracking is subject to errors in the location of the coordinates of the particles. Combined with polydispersity in the particle size distribution, this can influence structural measurements as we carry out here. In the case of 2-body correlation functions, the effect of polydispersity and tracking errors can be similar to a convolution with a Gaussian.88 In the case of amorphous systems, the effect of (mild) polydispersity has been investigated and this was found to have only a minor effect on the higher-order structure.89 In Fig. 10 in the appendix, we show the effects of adding a range of tracking errors to coordinates of computer simulation data. The height of the peaks and their width broadens as the error is increased.
![]() | (5) |
To reproduce the (small) polydispersity in the experimental system, we follow some previous work and use an equimolar five-component system.75 Here we set the diameters of the particles as σ1≤i≤5 = (0.93675445, 0.968377225, 1, 1.031622775, 1.06324555). We do not treat the effect of polydispersity on the dipolar interactions, other than rescaling the interaction range with the mean of the particle diameters. To our knowledge, no such study has been undertaken. We therefore implement the effects of polydispersity through the hard core contribution which we treat as follows.
To reproduce the (nearly) hard sphere behaviour of the experimental system, we use the Weeks–Chandler–Anderson (WCA)93 potential. Between particles i and j, takes the form:
![]() | (6) |
We added the dipole–dipole interaction shown in eqn (1), the Ewald sum for which is implemented with the KSpace package in LAMMPS.90 Here γ = γsim controls the strength of the dipolar interaction. The interaction potential for the simulations then reads
usim(r, θ) = uwca(r) + udip(r, θ). | (7) |
Here we quote simulation results in reduced Lennard-Jones units, that is to say the unit of length is the median diameter σ3, and time is in units of where m is the mass of a particle. We use the Barker–Henderson effective hard sphere diameter of the WCA component of the interaction to determine the effective volume fraction ϕeff in order to match the experiments, i.e.
where V is the volume, ni is the number of particles and σi,eff is the effective hard sphere diameter of the ith species. Each simulation run includes 1000 particles for ϕ = 0.1 and 3000 for ϕ = 0.3.
To prepare the system, the potential energy of random coordinates is minimised under eqn (7) to remove overlaps between particles. Then the system is evolved according to the Langevin dynamics. To match the timescales to the experimental system when it is out of equilibrium, we determine the Brownian time in the simulations τsimB = 0.979 ≈ 1. We therefore run the system for 200 time units prior to measurement which matches well to the 20 minutes in the case of the experiments. For equilibrium states (γ < 10, 12) for ϕ = 0.1 and 0.3 respectively, we equilibrate the simulations for 100 time units.
We also consider the 3-body spatial correlation function g3. Now this depends on the positions of three particles 1, 2, 3, i.e. three vectors, e.g. g3(r12, r23, r31) with the numbers reflecting the three particles. Here we elect to simplify our representation to the case where we fix two of the distances (to the particle diameter σ) so that the 3-body correlation function is plotted as a function of angle between them g3(r12 = σ, r23 = σ, η) where η is the angle between r12 and r23.
utcc(r, θ) = ulj(r) + udip(r, θ). | (8) |
Those resulting clusters which are rigid68 are shown in Fig. 1. In the case of zero field, we have the minimum energy Lennard-Jones clusters 8B, 9B, 10B, 11C and 12B.62,95 Upon application of the field, the clusters elongate in the z-direction. These stretched clusters are denoted S and are based on five-membered rings (9S, 10S, 11S and 12S). Further application of the field leads to clusters which are based on the 6A octahedron. These are denoted O (11O and 12O). Increasing the field still further leads to spiral clusters denoted P (9PAA and 10PAA). Additional minimum energy rigid clusters of particles interacting according to eqn (8) exist,68 but here we focus on those we find in our experiments and simulations. In our analysis of the TCC clusters, we set the Voronoi parameter fc = 0.82.67 Further details are of the identification of the clusters with basin hopping and inplementation in the TCC are available in ref. 68.
![]() | (9) |
Upon increasing the field strength such that γ ≈ 12 and ≈10 for volume fraction ϕ = 0.1 and 0.3 respectively, the string fluid becomes metastable to fluid-body-centered tetragonal phase coexistence.26 Since the system is now out of equilibrium, it is possible that small structural differences between the experiments and simulations may emerge, such as the small increase in 〈cos2θ〉 in the experiments with respect to the simulations for ϕ = 0.1.
As above in Section 4.1, at high field strengths when the system falls out of equilibrium, we see some discrepancy between experiment and simulation. In particular, we see stronger peaks, i.e. more ordering, in the simulation data for γ = 12 and 19 for ϕ = 0.1. We believe that this difference is too large to be attributed to tracking error. Interestingly, the difference between experiment and simulation is rather less significant in the case of ϕ = 0.3. We return to consider possible influences in the structure out of equilibrium in Section 6 below.
Since our system is anisotropic due to the external field, we expect to see some corresponding differences in structure between the plane perpendicular xy and the direction parallel z to the field. To explore this we now consider the pair correlation function in the xy plane g2xy(r) and z direction g2z(z). We expect that the formation of the string fluid may lead to significant structuring in the field direction, while in the perpendicular xy plane, there is repulsion between particles exactly in plane, but out-of-plane attractions lead to aggregation of the strings34 and ultimately the formation of the bct crystal. In Fig. 3(c) and (d), we show the pair correlations in the perpendicular plane for volume fraction ϕ = 0.1 and 0.3 respectively. In the case of ϕ = 0.1 and 0.3 at zero field strength we see good agreement between experiment and simulation. For a weak field at low volume fraction, there is some evidence of enhanced repulsion in the experiments with respect to the simulations, but this is not seen at higher volume fraction. Upon increasing the field strength such that γ = 35, we see a reasonable agreement between experiment and simulation. This is echoed in the higher volume fraction case (γ = 35) with again stronger peaks in the experiment. At the highest field strength, γ = 60 and 54 for ϕ = 0.1 and 0.3 respectively, the agreement between experiment and simulation is rather good.
Turning to the case of the correlations along the field direction, g2z(z) [Fig. 3(e) and (f)], for zero field strength, we recall that we may expect to see evidence of significant ordering upon application of the electric field, as the string fluid develops. At zero and low field strength, correlations are comparatively weak for both ϕ = 0.1 and 0.3 [though note the scale on the y-axis in Fig. 3(e) and (f)] and experiment and simulation are well matched. At higher field strength we find significantly stronger peaks in the case of the simulations for both volume fractions. This may be due to tracking errors, where we note the resolution of the confocal microscope is reduced in the z direction and thus tracking errors may be more severe here.
Summarising the behaviour we have uncovered through our analysis of the pair correlations, we find overall reasonably good agreement between experiment and simulation We suggest that discrepancies between experiment and simulation may be attributed to particle tracking errors (particularly in the case of g2z(z) at higher field strength).81,88
For both ϕ = 0.1 and 0.3, at zero and low field strength, a large and broad peak appears at η ≈ 60° which is consistent with an isotropic system where the interaction is angle independent.97 This is more prevalent in our experiments than in the simulations. Upon applying the field, a peak at η ≈ 180° emerges, relating presumably to the development of the string fluid. For ϕ = 0.1, γ = 12, This is more prevalent in the simulations than in the experiments, whereas at the slightly weaker field γ = 9 for ϕ = 0.3, the experiments exhibit a higher peak. When the system falls out of equilibrium, for both volume fractions, for we see rather enhanced ordering in the simulations for γ = 35 with a much stronger peak around η ≈ 120°. This may be related to a more rapid demixing towards a bct crystal. The peak towards η ≈ 180° is in fact weaker in the simulations than the experiments. For ϕ = 0.1, this is consistent with a slightly higher development of string-like structure in the experiments. At the highest field strength, both experiments and simulations show an increased degree of ordering, as evidenced in the peaks at η ≈ 60°, 120° and 180°. Again, the peaks seem somewhat more pronounced in the simulations.
To facilitate comparison between experiment and simulation, in each figure we fix the number of particles in the cluster in question. Now the topology of the minimum energy cluster changes upon increasing the dipolar contribution.68 In Fig. 6, we consider 8 and 9-membered clusters for volume fraction ϕ = 0.1 and 0.3. For larger clusters, there are fewer statistics at lower volume fraction and therefore, we focus on ϕ = 0.3 in Fig. 7. In line with the discussion in Section 3.6, and previous work,45,69,82 despite the lack of attraction between the hard spheres in the case for zero field strength, we nevertheless expect to find some clusters due to packing effects.
![]() | ||
Fig. 6 Populations of smaller minimum energy clusters detected by the topological cluster classification as a function of dipole strength. The number of particles detected in a given cluster Nc is scaled by the number of particles in a cluster of that size, Nm. Data points denote experiment and lines are computer simulation. Error bars in the experimental data represent the standard error. Shading denotes the change in cluster topology which minimizes the energy at different values of the dipole strength γ as indicated.68 White regions of graphs denote values of γ where the minimum energy clusters are not rigid, that is to say, single-particle or two-particle wide strings are formed, or other non-rigid structures.68 Data are shown for different cluster sizes m and volume fractions as follows. (a) Cluster size m = 8, volume fraction ϕ = 0.1. (b) m = 8, ϕ = 0.3. (c) m = 9, ϕ = 0.1. (d) m = 9, ϕ = 0.3. |
![]() | ||
Fig. 7 Populations of larger minimum energy clusters detected by the topological cluster classification as a function of dipole strength. As in Fig. 7, the number of particles detected in a given cluster Nc is scaled by the number of particles in a cluster of that size, Nm. Data points denote experiment and lines are computer simulation. Error bars in the experimental data represent the standard error. Shading denotes the change in cluster topology which minimizes the energy at different values of the dipole strength γ as indicated.68 White regions of graphs denote values of γ where the minimum energy clusters are not rigid. Data are shown for different cluster sizes m. (a) Cluster size m = 11, (b) m = 12. Here volume fraction ϕ = 0.3. |
Naively, one might expect that zero and low field strength would correspond to minimum energy clusters for the Lennard-Jones interaction as shown in Fig. 1 and known from studies with hard spheres,46,75,81,82 and that increasing the field strength might lead to a cascade of clusters of increasing elongation as indicated in Fig. 1. This turns out to be the case.
For the smallest size of cluster we consider, m = 8, indeed we see this trend with the Lennard-Jones minimum energy cluster 8B giving way to the 8O which minimises the potential energy for the dipolar system for both volume fractions [Fig. 6(a) and (b)]. Notably, experiment and simulation appear reasonably well-matched in the lower field case that the system is in equilibrium (γ ≤ 12, ϕ = 0.1 and γ ≤ 10, ϕ = 0.3). By this we mean that discrepancies are well within an order of magnitude (note that, for an 8-membered cluster to be successfully identified, failure to identify only one of the 8 particles will lead to the cluster not being identified, so “agreement” between experiment here is inevitably rather more stringent than in the case of pair correlations, say). The more significant discrepancy between the experiments and simulations emerges. At higher field strength, γ > 30 in both volume fractions ϕ = 0.1 and 0.3, we see an increase in the 8B modified pentagonal bipyramid population in the experiments compared to the simulations. Now this structure exhibits fivefold symmetry, and as such, is associated with non-crystalline ordering.61 We have noted above that, in the non-equilibrium conditions at higher field strength, the simulations appear more ordered. If the ordering in the simulations is crystalline, then a lack of five-fold symmetry with respect to the experiments would seem to be reasonable.
In the case of 9-particle clusters, the situation is more complex, as we consider three clusters, 9B, 9S and 9PAA. We find it instructive to start our analysis with ϕ = 0.3 [Fig. 6(d)]. At low field strength, γ = 0 and 9, there is a rather high population of 9S and 9PAA in both experiment and simulation. This dominates over the minimum energy cluster for the Lennard-Jones system, 9B which has two five-membered rings.67 At higher field strengths, the population of 9B drops rather precipitously in the simulations, but less so in the experiments. This is consistent with the case of m = 8 above, where the 8B cluster which, like the 9B, has a degree a fivefold symmetry is preferred. Notably though, the 9PAA dominates at all field strengths, which is not expected from energy considerations as it is the minimum energy cluster only for 12 ≲ γ ≲ 21. Now the 9PAA cluster is polytetrahedral in structure, and is enlongated with respect to the 9B. This is consistent with it being intermediate between the compact 9B in the case of zero field and strings of particles in the case of a strong field. From a structural perspective the 9S stretched polytetrahedron is intermediate between 9B and 9PAA in that it is the minimum energy structure for 7 ≲ γ ≲ 12. Its population is also intermediate between 9B and 9PAA for both experiments and simulations. For ϕ = 0.1 [Fig. 6(c)], there are rather fewer clusters identified, and none at all at low field strength. It is often the case that there are fewer (larger) clusters at lower volume fraction,70,72,81,98 so this in itself is reasonable. Like the case of ϕ = 0.3, the 9PAA dominates at all field strengths. At high field strength, we see more of the clusters with five-membered rings, 9B and 9S in the experiments than in the simulations.
We now consider larger clusters in Fig. 7(a)–(c). The statistics for larger clusters at the lower volume fraction are rather poor and therefore here we focus on a volume fraction ϕ = 0.3 respectively. Turning to 10-membered clusters [Fig. 7(a)], we see a somewhat similar behaviour to the m = 9 case. For both experiments and simulations, we see some of the Lennard-Jones minimum energy cluster, the defective icosahedron 10B at low field strength as one might expect. Like the 9-membered clusters, the 10B is present only in small quantities with the 10PAA dominant at all field strengths and 10S intermediate. In experiments the population of 10B (which exhibits some fivefold symmetry) is higher in the experiments than in the simulations, echoing the case for 8B in Fig. 6(a) and (b). Meanwhile the populations of 10S and 10PAA are comparable between experiments and simulations.
For 11-membered clusters, at low field strength, in experiment, we find 11S, 11SB and 11O dominating with rather less 11C, the latter being the minimum energy Lennard-Jones cluster, like the smaller Lennard-Jones clusters the 11C has a number of five-membered rings which may underly its small population. At higher field strength 11O dominates, similar to the case for 8O in m = 8. However for m = 11 the additional two structures 11S and 11SB seem to lie between 11C and 11O, and their population falls away at higher field strength. This is broadly consistent with expectations of 11O dominating at higher field strength (it is the minimum energy cluster for 27 < γ < 38).
Finally, the 12-membered clusters are shown in Fig. 7(b). Here in the experimental data, 12O, which is the minimum energy cluster for 19 < γ < 31 is the most popular cluster at all field strengths. We see some 12B and 12SB at weak field strengths, but these vanish for field strengths greater than γ = 15. The simulations exhibit the same qualitative trend, of 12O dominating.
We therefore probe the orientation of some dipolar-Lennard-Jones clusters using our method described above in Section 3.7. The principle is illustrated in Fig. 8(a) and (b). Here, renderings of experimental data where clusters identified as 10PAA for volume fraction ϕ = 0.3 are shown. In the case of zero field strength [Fig. 8(a)], no preference in cluster orientation is seen. For γ = 54 [Fig. 8(b)], we see that the clusters have oriented with the field. We now explore this phenomenon quantitatively. In Fig. 8(c) and (d), for volume fractions ϕ = 0.1 and 0.3 we plot the degree of alignment with the field for experimental data (data points) and simulation (lines).
![]() | ||
Fig. 8 Orientation of anisotropic clusters with the electric field. (a) shows a 3d plot of the clusters (10PAA) found in the experiments at ϕ = 0.3 at zero field, where the principal axis (indicated as the black line) of each cluster does not align with the electric field. Whereas in (b) these 10PAA clusters showed a higher degree of alignment along the field (z-axis) as the field takes a value of γ = 40. (c) and (d) 〈P2(cosψ)〉 of the principal axis is plotted as a function of reduced field strength, for the anisotropic clusters 8O, 9PAA, 9S, 10PAA, 11O and 12O. for the two volume fractions under consideration. (c) ϕ = 0.1 and (d) ϕ = 0.3. Data points are experimental data and lines are computer simulation. Shading denotes the system falling out of equilibrium.26 |
We begin our analysis with ϕ = 0.1 [Fig. 8(c)]. As exemplified by 8O (for which we have the best statistics), there is a trend towards more alignment with the field as its strength is increased. The experimental data in general shows an increase in alignment, while for the simulations, 9PAA, 10PAA and 12O show some non-monotonicity in their response to the field. Turning to the case of higher volume fraction, in Fig. 8(d), we see a similar behaviour. Experimental data all show an increase in alignment with the field. Contrasting with this, in the simulations, 11O even seems to approach to 〈P2(cosψ)〉 = −1/2, the case for alignment in the xy plane. This may be due to the cluster being identified in the bct crystal in which it may lie perpendicular to the field.
![]() | ||
Fig. 9 Aging in the phase separation regime. Here we consider the state point ϕ = 0.3, γ = 32. (a) Radial distribution functions g2(r) are plotted for the times indicated. These are compared with experimental data (for a waiting time around 200τB). Experimental data [plotted in Fig. 3(a)] and simulations data taken for t = 0, 30, 100, 300, 900 time units. (b) Time-evolution of selected TCC clusters for the simulation. Here time is in Lennard-Jones time units. |
We now turn to the higher-order correlations and in particular plot the time-evolution of certain TCC clusters as a function of time in Fig. 9(b). These show rather clear trends following their geometry. We have remarked that the clusters incorporating octahedra, in particular 8O and 12O are compatible with the bct crystal and these are indeed found at high field strength [Fig. 5(b)]. The population of these clusters increases with time here which is consistent with the system moving towards the fluid-bct state. Conversely the Lennard-Jones cluster 8B and 9S which also exhibits some five-fold symmetry find their population falling, presumably because their symmetry is incompatible with the fluid-bct state. Finally, clusters which are ground states at intermediate field strength, 9PAA and 10PAA68 exhibit relatively little change in their population. This suggests that the TCC may be a suitable means to probe the time-evolution of the dipolar colloids. A more extensive investigation we leave for the future.
(i) Fig. 2 shows confocal microscopy images of our system at volume fraction ϕ = 0.1, taken along the xy plane (perpendicular to the direction of the electric field) and yz plane (along the direction of field) for both zero dipole strength and at maximum dipole strength at γ = 45. We see string formation along the direction of the field at γ = 45.
Now our system becomes metastable to fluid-bct crystal phase coexistence at γ ≥ 10 and 12 for volume fraction ϕ = 0.1 and ϕ = 0.3 respectively. Under these conditions, given that throughout most of this work, we do not treat the time-evolution of the system, discrepancies between the experiments and simulations are possible due to each taking a different route through the energy landscape. However, our analysis in Fig. 9 suggests that for two-point correlation functions, the effect of the time-evolution is quite weak. For the TCC clusters, we do see some changes, and we discuss this in more detail below.
Fig. 2(c) shows a plot of bond order parameters of the string fluid from both our experiment and simulation results, similar to the study published by Li et al.30 Here, our simulations (line) and experiments (data points) show good agreement. As the field strength increases, the bond angle θ tends towards 180°. As indicated in Fig. 2(c), the degree of string formation increases which is consistent with the string formation is continuous rather than a sharp transition.
(ii) In Fig. 3 we plot g2(r). This is generally in reasonable agreement between computer simulations and experiments for both ϕ = 0.1 and ϕ = 0.3. We can therefore be fairly confident that the simulation model used in our work is a reasonable reflection of our experimental system. We see the emergence of long range order as field strength is increased. This is expected since the confocal images show that as the fluid becomes more structured as the colloids aligned along the field when it is switched on. The pair correlation function g2 does however show significant discrepancies emerge at high field strength. The increased ordering in the simulations may be due to the system being further down the path to forming a bct crystal than is the case in experiment.
The results of g2xy(r) and g2z(z) are both consistent with that of g2(r) where the height of the first peak increases as field strength is switched on. The significant differences in peak positions and shapes between g2xy(r) and g2z(z) show the fluid structure across the xy-plane differs from fluid structure along the z-axis. We also observe peak splitting (g2xy(r)) from a broad peak into two distinct peaks as field strength increases from γ = 35 to γ = 60/54, showing an increase in ordering of fluids across the xy-plane. Whereas ordering in the fluid along the z-axis occurs above γ = 12/9. However, fluid structures still show little difference at low field strength, except that the first peak in particular is stronger in the case of the simulations, which we attribute to particle tracking errors. At higher field strengths as the system falls out of equilibrium and starts to order, larger discrepancies emerge, notably in the higher first peak in g2xy(r) and much higher peaks in gz(z). The latter we attribute to tracking errors.
(iii) We now consider three-body correlations in Fig. 4. As the dipolar strength increases, the peak at 180° increases as expected since more dipolar colloids form strings. However, we also observe peaks at 60° and 120° increasing with respect to field strength. The triplet correlations in the simulations show increasing structure at higher field strength with respect to the experiments.
(iv) The renderings in Fig. 5 and the plots in Fig. 6 and 7 show the population of dipolar clusters of different geometries and sizes analysed with the TCC. Overall, we find that the clusters we observe in our system follow reasonably those of the dipolar-Lennard-Jones clusters (Fig. 1). That is to say, we see more elongated clusters at higher field strengths. In all cases, at high field strength it is the LJ-dipolar cluster that corresponds to the highest dipolar interaction that we find in both experiment and simulation. We find this to be a significant outcome of this work, providing strong evidence in support of modelling colloids in an AC electric field with dipolar interactions.
Minimum energy dipolar clusters by definition imply zero temperature, and are determined by the interaction energy. However, at finite temperature, entropy plays a role. The fact that the cluster population trends largely follow the minimum energy clusters indicates that the behavior of dipolar colloids is strongly influenced by energetics. This is in stark – and surprising – contrast to earlier work which suggested that energy plays only a very limited role in observed cluster populations.71 That work investigated Frank's well-known conjecture that icosahedra “will be a very common grouping in liquids”.61 In fact, at the triple point of the Lennard-Jones system, only one particle in 1000 was found to be in an icosahedron and other 13-membered clusters dominate,71 quite unlike the findings here in which the minimum energy structure dominates at high field strengths. Presumably the strength of the dipolar interactions (which are much larger than e.g. interactions in the Lennard-Jones system when it is in the liquid state71) is important here. At high field strength, for m = 8, 9, 10 clusters we find more clusters with five-membered rings in the experiments than in the simulations. This we attribute to geometric frustration as the system falls out of equilibrium, with experimental and simulated realisations of the system taking different paths in the energy landscape (see below).
(v) Fig. 8 shows that some anisotropic dipolar clusters tend to align along the z-axis (parallel with direction of the field) while others in fact orient perpendicular to the field. This is a surprising result and we attribute it to incipient crystallisation in the system. Inspection of Fig. 5 shows ordering and a structure dominated by 8O and 12O. While in the experiments the 8O and 12O rich regions are oriented vertically, it is possible that in simulation a somewhat different path through the energy landscape is followed.
(vi) Fig. 9 provides an indication of how the time-evolution of the system may be monitored. We see that changes in the radial distribution function are limited, but the time-evolution of the TCC clusters is very much more marked. This suggests that this sort of analysis may be useful to probe kinetic pathways in the dipolar colloidal system.
We now wish to discuss the relevance of our results with higher order structure. Higher-order structure is sensitive to changes in structure in a way the two-point correlation function are not. For example, at low volume fraction, even relatively weak field strengths (γ = 12) cause a dramatic change in the population of, for example the 8O cluster [Fig. 6(a)] while g2(r) exhibits rather little change. We can therefore conclude that higher order analysis such as that presented here is much more sensitive to small variation in structure and interactions. Comparing our study with other previously published work which used a TCC analysis of gels and glasses,38,89,100–103 we can conclude that such higher order structural analysis is better at capturing the onset of structural changes in amorphous systems than pair correlations g2(r), as may be inferred from other work.37,39,104,105 What is new here is that we have considered a system with anisotropic interactions.
The more significant discrepancies that we find are in the regime in which the system departs from equilibrium and becomes metastable to fluid-bct crystal phase coexistence. Now the early stages of this transition have been investigated recently.34 However, in related phenomena, such as the condensation of colloids with an effective attraction to form a gel network, the role of hydrodynamic interactions was found to be very important. Hydrodynamics control not only the timescale for the condensation,106 but also the higher-order structure of the resulting non-equilibrium gel network.48,107 In (non-equilibrium) gelation of particles with spherically symmetric attractions, experiments exhibit many fewer clusters with fivefold symmetry than do Brownian dynamics computer simulations,48 quite the opposite trend of what is observed here in Fig. 6 and 7. We suggest that careful study, using simulations with hydrodynamic interactions and time-resolved experimental observation along the lines of ref. 48 may enable an analysis of the time-evolution of this system more closely matched to the experiments than we have been able to perform here. This would be very interesting to probe in the future.
![]() | ||
Fig. 10 The effects of adding errors to the coordinates on the radial distribution function g(r). Shown is the data for a volume fraction ϕ = 0.1 and field strength γeff = 50. Experimental data (circles) and data without errors added to the simulated coordinates (red line) are reproduced from Fig. 3(a). Errors are added to the simulated coordinates using a Gaussian with standard deviation 0.05 (pink line) 0.1 (grey line) and 0.2 (black line). |
Footnote |
† A value of 3.8% was obtained for the polydispersity using static light scattering, and a higher value was obtained using electron microscopy. Different methods are known to provide different results when measuring the size polydispersity.80,81 In this case, that the particles crystallize readily suggests that their polydispersity was around 5% or less. |
This journal is © The Royal Society of Chemistry 2025 |