Tuning higher order structure in colloidal fluids

Abstract

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.

---

1 Introduction

Particles with a dipolar interaction are of significant fundamental importance in the study of fluids and disordered materials. They are among the simplest models which describe long range directional interactions, which are exhibited by molecules.^1,2 Colloidal dispersions provide suitable models of atomistic and molecular systems as they exhibit phase behaviour following the same rules of statistical mechanics, yet are amenable to real space observation using optical microscopy.^3–5

Colloidal dipolar systems fall broadly into two categories. Some, for example ferromagnetic nanoparticle systems, like atoms and molecules, have an intrinsic dipole moment⁶ and can be modeled with the Stockmayer model which combines a dipolar interaction with a Lennard-Jones interaction.⁷ These systems exhibit intriguing string-like structures,^8,9 with branching,¹⁰ coiling¹¹ and clustering⁷ behavior, not to mention a ferromagnetic transition.¹² 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.¹⁴ 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 interactions^16,17 and direct observation of cluster growth.¹⁸ In addition to their fundamental interest, such dipolar colloidal systems may find application as electrorheological fluids,^19,20 hydraulic valves²¹ and photonic materials.²²

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 real²³ and reciprocal space^24,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.²⁰ Tuning the electric (or magnetic) field in situ enables the control of phenomena such as a martensitic transition.²⁷ Adding softness²⁸ or attractions²⁹ 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 analytically³¹ 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,¹⁸ and aggregation phenomena between the strings.³⁴ This suggests that the structure of the “string fluid” may be rather interesting and it has been investigated analytically³¹ and in reciprocal space⁹ 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 studies⁴ 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 arrest^35,37–39 and polymorph selection^40–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 g₃,^44,47,48 and higher-order correlations such as common neighbor analysis (CNA)⁴⁹ and Voronoi face analysis.⁵⁰ 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.,⁵³ 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 liquid^54,55 and also in the characterization of fivefold symmetric order in amorphous systems.³⁷ 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 supervised⁵⁹ and unsupervised learning⁶⁰ 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,⁶¹ 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,⁶² 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,⁶³ “Sticky Sphere”,⁶⁴ Stockmayer⁶⁵ models and more exotic interactions relating to multiple fields.⁶⁶ Since the dipolar interaction of the Stockmayer model is not constrained to lie in any particular direction it thus corresponds to a molecular (or nanoparticle⁷) 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,³⁸ colloid–polymer mixtures interacting via the Morse potential,^69,70 colloidal suspensions with attractive interactions,⁷¹ colloidal gels^48,72 and the liquid–gas interface.⁷³ 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),⁶⁸ 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.²⁶ We consider pair correlations in the form of radial distribution functions g₂(r) and three-body correlations in the form of order parameters to determine “string-like” configurations and also the triplet correlation function g₃(r,r′,η). We use the topological cluster classification^67,68 to explore higher-order spatial correlations in the form of minimum energy clusters of the dipolar-Lennard-Jones interaction.⁶⁸ Before concluding, we provide an outlook on the study of non-equilibrium behaviour accessed via higher-order correlations.

2 Dipolar interactions in colloidal systems

The colloids in our experiments are suspended in an index-matching solvent with added salt, as described in Section 3.1. Like very many colloidal systems,⁷⁴ ours carry an electrostatic charge. This we screen by adding salt, and we have previously demonstrated that this is sufficient to treat the particles as hard spheres, with a slight increase in effective diameter due to the residual electrostatics.⁷⁵

When colloids are subjected to an external electric field E, a dipole–dipole interaction is induced between the particles which reads

where u_dip is the dipolar interaction, k_BT is the thermal energy and θ is the angle made by r and the z-axis. Here r is the vector between the centres of the colloids. In our experimental system, γ = γ_exp is a dimensionless prefactor that depends on the strength of the external field and material properties of the system. Here r is the vector connecting the centres of the two particles. The prefactorwhere ε_s is the dielectric constant of the solvent and the dipolar momenthere α = (ε_p − ε_s)/(ε_p + 2ε_s) where ε_p and ε_s is the dielectric constant of the particles and solvent respectively and E_loc is the local electric field. Here we follow Hynninen and Dijkstra²⁶ and set E_loc = E/(1 − απ/6) which is the case for a cubic crystal.

Combining with the hard sphere interaction noted above u_hs, the total interaction between two colloids under the electric field becomes

u_total(r, θ) = u_hs(r) + u_dip(r, θ). (4)

3 Methods

3.1 Experimental

The colloidal suspension used in this experiment was prepared by adding sterically-stabilised polymethyl methacrylate (PMMA) spheres (synthesized following ref. 76 and 77) (ρ = 1.196 g cm⁻³)⁷⁸ of diameter σ = 1.73 μm and polydispersity ≲5%⁷⁹ † to a mixture of density and refractive index matched solvents. The particles we labelled with 1,1′-dioctadecyl-3,3,3′,3′-tetramethyl-indocarbocyanine perchlorate, which may be excited at a wavelength of 543 nm.⁷⁶ The solvent is a mixture of cis-decalin (ρ ≈ 0.897 g cm⁻³) and cyclohexyl-Bromide (CHB) (ρ = 1.32 g cm⁻³). Tetrabutylammonium bromide (TBAB) salt was dissolved in the solvent to make up a solution with TBAB concentration of 260 μM. This corresponds to a Debye length κ⁻¹ of around 100 nm.⁸² Since the Debye length is much less than the particle diameter, in the absence of an electric field, the colloids behave as nearly hard spheres.^74,82 While more sophisticated treatments may be carried out to match the interaction potential,^74,82 here we use a slightly soft potential in the computer simulations and presume this to be sufficient to match the experimental system, noting that the effects we seek to study are dominated by the dipolar interactions, rather than the hard core or precise value of the colloid volume fraction (which we determine by weighing out the samples).⁸⁰

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 ⁸³ 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.⁴ 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⁻¹ measured peak-to-peak. Before each measurement, the system is allowed to stabilise for at least 20 min. The Brownian time τ_B = (3πνσ³)/(4k_BT) ≈ 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.

3.2 Particle tracking

Here we use a slight modification to enhance the accuracy of the coordinates that we detect.⁸⁴ We begin by carrying out a conventional centroid location.^85,86 This seeks the brightest pixels and weights the brightnesses of the surrounding pixels to obtain an estimate of the centre of the colloidal particle. Overlaps corresponding to multiple pixels within a single particle being identified as separate particles are removed. This method works well in (quasi) 2d studies,⁸¹ but in the case of the 3d confocal microscopy that we carry out here, overlaps between blurred images of particles in the z-direction can be a problem.

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.⁸⁷

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.⁸⁸ 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.⁸⁹ 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.

3.3 Computer simulation

We employed computer simulations using the LAMMPS package in the NVT ensemble.⁹⁰ To mimick the colloidal motion, we used Langevin dynamics, where the ith particle has position r_i and velocity v_i and the velocities evolve in time aswhere m is the mass of the particle, V is the total potential energy while μ is a friction constant and ξ_i a random noise force. The friction constant sets the time scale for the decay of velocity correlations as t_d = m/μ. For colloids, the physical situation corresponds to an overdamped limit with t_d very small. It has been found that setting t_d = 0.1 here is appropriate to following the dynamics of non-equilibrium colloidal systems and can robustly be compared with experiments.^91,92

To reproduce the (small) polydispersity in the experimental system, we follow some previous work and use an equimolar five-component system.⁷⁵ 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)⁹³ potential. Between particles i and j, takes the form:

where ε_wca = 10k_BT is the interaction energy and σ_ij = (σ_i + σ_j)/2.

We added the dipole–dipole interaction shown in eqn (1), the Ewald sum for which is implemented with the KSpace package in LAMMPS.⁹⁰ Here γ = γ_sim controls the strength of the dipolar interaction. The interaction potential for the simulations then reads

u_sim(r, θ) = u_wca(r) + u_dip(r, θ). (7)

Here we quote simulation results in reduced Lennard-Jones units, that is to say the unit of length is the median diameter σ₃, 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, n_i 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 τ^sim_B = 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.

3.4 Bond order parameters for dipolar colloids

One method to quantify the angular correlation between particles as a function of external electric field strength is the string fluid order parameter.³⁰ This has already been shown to be sensitive to variation in field strength by Li et al.³⁰ This is calculated by finding the angle θ made by a reference particle with its two nearest neighbors [see Fig. 2(c) inset]. Here we take 〈cos² θ〉 as the order parameter for the string fluid. For a system consisting of perfect strings, 〈cos² θ〉 = 1.

Fig. 2 The colloidal dipolar system. (a)–(d) Representative images using confocal microscopy in the horizontal xy (a) and (b) and vertical xz (c) and (d). Here volume fraction ϕ = 0.1. (a) and (c) No field, E = 0. (b) and (d) Confocal microscopy images of a system at the same volume fraction at the maximum electric field strength corresponding to γ = 60. (e) String fluid order parameter (〈cos² θ〉) as a function of external electric field strength E. The inset indicates the angle θ. Data are shown for experiments (data points) and simulations (lines) for volume fractions ϕ = 0.1 and 0.3 as indicated. Scale bars in (a)–(d) are 20 μm. Error bars are the standard error.

3.5 Two- and three-body correlation functions

We also calculated the two-body spatial correlation function, the radial distribution function g₂(r). In addition to the isotropic g₂(r), we consider g_2xy(r) which measures correlations in the xy plane perpendicular to the field and g_2z(z), which measures correlations in the z direction along the field. For g_2xy(r), only pairs of particles that are perpendicular to the xy plane with a tolerance of ±5° were used. In the case of g_2z(z), particle pairs that are parallel to the z-axis were chosen with a ±5° tolerance.

We also consider the 3-body spatial correlation function g₃. Now this depends on the positions of three particles 1, 2, 3, i.e. three vectors, e.g. g₃(r₁₂, r₂₃, r₃₁) 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 g₃(r₁₂ = σ, r₂₃ = σ, η) where η is the angle between r₁₂ and r₂₃.

3.6 Topological cluster classification

As discussed in the introduction, the topological cluster classification identifies local geometric motifs whose bond topology (defined here through a modified Voronoi decomposition) is identical to that of minimum energy clusters of a specific size.⁶⁷ These clusters are identified using the energy optimization algorithm GMIN which uses basin-hopping to find the local energy minimum that corresponds to a specific configuration for a number of particles in isolation.⁶³ Now such minimum energy clusters require an attractive interaction, and therefore to investigate the effect of the dipolar interaction, clusters were determined for a dipole added to a Lennard-Jones interaction.⁶⁸ That is to say, the interaction potential for which the minimum energy clusters were determined was

u_tcc(r, θ) = u_lj(r) + u_dip(r, θ). (8)

where u_lj is the Lennard-Jones (LJ) interaction. Although the attractive Lennard-Jones contribution is not part of the experimental (or simulated) system we consider here, the importance of packing effects on the structure of liquids and dense fluids has a long history.⁹⁴ The seminal work of Weeks, Chandler and Anderson⁹³ established the role of packing in the structure of dense fluids in the absence of attractions. It has recently been shown that the higher-order structure of dense hard spheres is closely related to those of attractive systems.^45,46 In fact, the higher-order structure of the Lennard-Jones and WCA systems, along with hard spheres, as determined through the topological cluster classification are rather similar.⁷¹ We therefore expect that the use of dipolar-LJ clusters will likewise be reasonable here and in any case will provide a suitable measure of the change in the fluid structure under the electric field.

Those resulting clusters which are rigid⁶⁸ 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,⁶⁸ 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 f_c = 0.82.⁶⁷ Further details are of the identification of the clusters with basin hopping and inplementation in the TCC are available in ref. 68.

3.7 Orientation of anisotropic clusters

Since the dipolar interaction is anisotropic, the clusters found by the TCC may have a preferred orientation with respect to the direction of the electric field. We consider clusters found in both experiments and simulations. We calculated the principal axis by taking the eigenvector of the inertia tensor of the cluster with the largest eigenvalue. We then determine the angle made by the principal axis with respect to the direction of the electric field. To quantify this angle distribution, we use an order parameter commonly used for liquid crystals,where ψ is the angle made by the principal axis of each cluster with the direction of the electric field.⁹⁶ In the isotropic case when the field is switched off, there should be no preferred orientation. In the perfectly aligned state, all the clusters should lie parallel to the electric field so this order parameter is 1. We note that alignment with the electric field assumes that the clusters are isolated and free to rotate. The system can exhibit a fluid-body centred tetragonal crystal phase coexistence. Some clusters (8O, 11O and 12O) are in fact compatible with this crystal structure and may lie perpendicular to the field, in which case P₂(ψ = π/2) = −1/2.

4 Results

We now present both experimental and simulation results for colloidal dipolar fluids at volume fractions ϕ = 0.1 and ϕ = 0.3. We investigate the bond order parameters, two and three body correlation functions, and populations of minimum energy clusters and the orientations of these clusters. These quantities are considered as a function of dipole strength γ. We convert external field strength to the dipole strength γ using eqn (2), which we use across experiments (eqn (3)), simulations (eqn (7)) and minimum energy clusters (eqn (8)).

4.1 Bond order parameter analysis of the string fluid

In our study, with the string fluid order parameter 〈cos² θ〉 (see Section 3.4), we explore slightly higher colloid volumes fraction (0.1 and 0.3) than some previous work.³⁰ Our results are shown in Fig. 2 and for the string order parameter, we find good agreement between experiment and simulation if we scale the values of γ in the experiment by 1.2, which we do henceforth. The value of 〈cos² θ〉 is small for both experiment and simulation in the case of zero field, and then shows a significant increase with field strength to a value of 〈cos² θ〉 ∼ 0.8 for ϕ = 0.1 and ∼0.65 for ϕ = 0.3.

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.²⁶ 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 〈cos² θ〉 in the experiments with respect to the simulations for ϕ = 0.1.

4.2 Pair correlation function

We continue our analysis by considering pair correlations in the form of the radial distribution function g₂(r) in Fig. 3. The orientationally averaged g₂(r) is shown in Fig. 3(a) and (b) for volume fraction ϕ = 0.1 and 0.3 respectively for both experiment (data points) and simulation (lines). At zero field strength, in the hard sphere limit, we see reasonable agreement between experiment and simulation (as has been noted previously).⁸¹ The slightly higher first peaks of the simulation data may be attributed to the particle tracking errors in the experiments.⁸⁸ For weak field strengths (γ = 9), again we see comparable agreement between experiment and simulation to that of the hard sphere case.

Fig. 3 Pair correlations in the colloidal dipolar system. Data are shown for both experiment (points) and computer simulations (lines) at volume fraction (a), (c) and (e) ϕ = 0.1 and (b), (d) and (f) ϕ = 0.3. Data are taken for different field strengths expressed through the parameter γ as indicated. Data are offset for clarity. (a) and (b) show g₂(r), (c) and (d) show g_2xy(r) by considering correlations in the xy plane, (e) and (f) show g_2z(r) where correlations are taken along the z axis.

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 g_2xy(r) and z direction g_2z(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 strings³⁴ 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, g_2z(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 g_2z(z) at higher field strength).^81,88

4.3 Three-body correlations

As noted above, particle-resolved studies lends itself to analysis of higher-order correlations^{35,37–39,81} and here we begin with the three-body correlation function g₃. This may be represented in a variety of ways and here we chose to show the dependence of g₃ upon the angle η between two particles with respect to a particle of interest as shown in Fig. 4(a). We set the distance between the two particles and the particle of interest to be the diameter, so that the quantity plotted is g₃(η).

Fig. 4 The three-body correlation function g₃(η). (a) Schematic indications of geometries of interest, with values of the bond angle η = 60°, 120° and 180° (b) g₃(η) is plotted for volume fraction ϕ = 0.1. (c) g₃(η) for volume fraction ϕ = 0.1. In (b) and (c), data points are experimental data and lines are simulation data.

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.⁹⁷ 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.

4.4 Populations of minimum energy clusters

We now move to still higher-order spatial correlations and consider the dipolar-Lennard-Jones minimum energy clusters identified by the topological cluster classification (Fig. 1).⁶⁸ In Fig. 5, to give an overview of the response of the cluster population to the electric field, we render two example state points. These are ϕ = 0.3, γ = 0 in Fig. 5(a) and ϕ = 0.3, γ = 54. We see some Lennard-Jones clusters in the case of zero field (8B, 9B and 11C). For the higher field strength, we see an organisation into columns (related to the gradual demixing towards the bct crystal). These are dominated by 8O and 12O clusters. We now proceed to analyse the cluster population in our system in more detail.

Fig. 5 Rendering of TCC clusters using experimental coordinate data. Key indicates the particular cluster a particle is in. Particles are rendered in the largest cluster in which they are identified. Small grey particles are not identified in any of the clusters considered. (a) ϕ = 0.3, γ = 0. (b) ϕ = 0.3, γ = 54.

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.⁶⁸ 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 N_c is scaled by the number of particles in a cluster of that size, N_m. 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.⁶⁸ 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.⁶⁸ 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 N_c is scaled by the number of particles in a cluster of that size, N_m. 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.⁶⁸ 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.⁶¹ 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.⁶⁷ 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.

4.5 Cluster orientation

We have observed that even at zero dipole strength, some of the minimum energy dipolar clusters (such as 10PAA at ϕ = 0.3) are present in our system [Fig. 6(f)]. Now the dipolar-Lennard-Jones clusters that we consider (Fig. 1) are aligned with the dipolar interactions, and thus with the electric field. It is therefore reasonable to suppose that the clusters might exhibit some alignment with the field, and that this would increase as a function of field strength.

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) 〈P₂(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.²⁶

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 〈P₂(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.

5 Time-evolution

We have remarked that in the regime that γ > 10, 12 for ϕ = 0.1, 0.3, the system becomes metastable to a phase separated fluid-bct crystal state.²⁶ As an outlook, we briefly consider the consequences of time-evolution for our analysis. To this end, we run a longer simulation of 1000 time units that we analysed from the beginning (without applying any equilibration time after the initial energy minimisation). We consider the ϕ = 0.3, γ = 32 state point. In Fig. 9(a), we plot the radial distribution function, and also reproduce the experimental data for this state point from Fig. 3(b). The time-evolution of the g₂(r, t) does show some changes, notably a split second peak can be detected longer times. This is compatible with increased ordering.⁹⁹ The changes in the pair-correlations are, however, rather small.

Fig. 9 Aging in the phase separation regime. Here we consider the state point ϕ = 0.3, γ = 32. (a) Radial distribution functions g₂(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 10PAA⁶⁸ 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.

6 Discussion

We now discuss our findings figure by figure.

(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.³⁰ 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 g₂(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 g₂ 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 g_2xy(r) and g_2z(z) are both consistent with that of g₂(r) where the height of the first peak increases as field strength is switched on. The significant differences in peak positions and shapes between g_2xy(r) and g_2z(z) show the fluid structure across the xy-plane differs from fluid structure along the z-axis. We also observe peak splitting (g_2xy(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 g_2xy(r) and much higher peaks in g_z(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.⁷¹ That work investigated Frank's well-known conjecture that icosahedra “will be a very common grouping in liquids”.⁶¹ 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,⁷¹ 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 state⁷¹) 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 g₂(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 g₂(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.³⁴ 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,¹⁰⁶ 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,⁴⁸ 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.

7 Conclusion

We have performed a detailed analysis of the string fluid structure in an anisotropic system of dipolar colloids and found reasonable agreement between our experiments and computer simulation data across a wide range of interactions tuned with the electric field. We found both bond-order parameter analysis of strings and the three-body correlation function g₃ to be suitable to quantify the degree of string formation in dipolar colloids but with g₃ offering more detailed information and can be used as a form of “colloidal finger-print”. Using the topological cluster classification, we find that our experiments and simulations broadly agree with expectations from minimum energy clusters of a dipolar-Lennard-Jones system.⁶⁸ That is to say, structural transformations predicted at zero temperature for a Lennard-Jones-Dipolar system are rather effective in their prediction of higher-order structure in the nearly-hard sphere-dipolar experiments and simulations. At high field strength, the cluster population in both our experiments and simulations is dominated by the minimum energy clusters for all sizes 8 ≤ m ≤ 12. Finally, not only can we identify clusters relevant to the dipolar system but also to investigate their orientation with respect to field strength.

Data availability

Data for this article, including image data and codes are available at Zenodo at https://doi.org/10.5281/zenodo.12788211.

Conflicts of interest

There are no conflicts to declare.

Appendix

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).

Acknowledgements

XW and FCM thank European Research Council (ERC) for funding an Advanced Grant for the DYNAMIN project, grant agreement No. 788968. CPR acknowledges the Agence Nationale de Recherche for grant DiViNew. We thank the KITP for hospitality where some of this work was done. YY gratefully acknowledges the China Scholarship Council. FJM was supported by a studentship provided by the Bristol Centre for Functional Nanomaterials (EPSRC Grant No. EP/L016648/1). The authors would like to thank Mark Miller, Josh Robinson and Peter Crowther for the help in building the TCC software. The authors would also like to thank Mike Allen and Didi Derks for helpful discussions.

References

1. J. N. Israelachvili, Intermolecular and Surface Forces, Academic Press, 3rd edn, 2011.
2. M. A. Gebbie, A. M. Smith, H. A. Dobbs, A. A. Lee, G. G. Warr, X. Banquy, M. Valtiner, M. W. Rutland, J. N. Israelachvili, S. Perkin and R. Atkin, Chem. Commun., 2017, 54, 1214–1224.
3. R. Evans, D. Frenkel and M. Dijkstra, Phys. Today, 2019, 72, 38.
4. A. Ivlev, H. Löwen, G. E. Morfill and C. P. Royall, Complex Plasmas and Colloidal Dispersions: Particle-resolved Studies of Classical Liquids and Solids, World Scientific Publishing Co., Scientific, Singapore, 2012.
5. B. Bharti and O. D. Velev, Langmuir, 2015, 31, 7897.
6. M. A. Boles, M. Engel and D. V. Talapin, Chem. Rev., 2016, 116, 11220–11289.
7. E. V. Novak, E. S. Pyanzina, P. A. Sánchez and S. S. Kantorovich, Eur. Phys. J. E: Soft Matter Biol. Phys., 2019, 42, 158.
8. M. Klokkenburg, R. P. A. Dullens, W. K. Kegel, B. Erné and A. P. Philipse, Phys. Rev. Lett., 2006, 96, 037203.
9. E. Pyanzina, S. S. Kantorovich, J. J. Cerdà, A. O. Ivanov and C. Holm, Mol. Phys., 2009, 107, 571–590.
10. S. S. Kantorovich, A. O. Ivanov, L. Rovigatti, J. M. Tavares and F. Sciortino, Phys. Chem. Chem. Phys., 2015, 17, 16601–16608.
11. A. Spatafora-Salazar, S. Kuei, L. H. P. Cunha and S. L. Biswal, Soft Matter, 2023, 19, 2385–2396.
12. J. Weis and D. Levesque, J. Chem. Phys., 2006, 125, 034504.
13. J. G. Donaldson, P. Linse and S. S. Kantorovich, Nanoscale, 2017, 9, 6448–6462.
14. D. Mostarac, E. V. Novak and S. S. Kantorovich, Phys. Rev. E, 2023, 108, 054601.
15. N. Elsner, D. R. E. Snoswell, C. P. Royall and B. V. Vincent, J. Chem. Phys., 2009, 130, 154901.
16. E. V. Yakovlev, K. A. Komarov, K. I. Zaytsev, N. P. Kryuchkov, K. I. Koshelev, A. K. Zotov, D. A. Shelestov, V. L. Tolstoguzov, V. N. Kurlov, A. V. Ivlev and S. O. Yurchenko, Sci. Rep., 2017, 7, 13727.
17. V. Yakovlev, N. P. Kryuchkov, S. A. Korsakova, N. A. Dmitryuk, P. V. Ovcharov, M. M. Andronic, I. A. Rodionov, A. V. Sapelkin and S. O. Yurchenko, J. Colloid Interface Sci., 2021, 608, 564–574.
18. A. Ben Salah, T. Ukai, L. Mingyuan, H. Morimoto and T. Maekawa, AIP Adv., 2023, 13, 055010.
19. W. M. Winslow, J. Appl. Phys., 1949, 20, 1137–1140.
20. U. Dassanayake, S. Fraden and A. van Blaaderen, J. Chem. Phys., 2000, 112, 3851.
21. M. Parthasarathy and D. Klingenberg, Mater. Sci. Eng., 1996, R17, 57.
22. A. van Blaaderen, MRS Bull., 2004, 29, 85–90.
23. A. Yethiraj and A. van Blaaderen, Nature, 2003, 421, 513.
24. M. Klokkenburg, B. Erné, A. Weidman, A. V. Petukhov and A. P. Philipse, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2007, 75, 051408.
25. A. A. Wiedenmann, U. U. Keiderling, M. Meissner, D. D. Wallacher, R. Gähler, R. P. May, S. Prévost, M. Klokkenburg, B. Erné and J. Kohlbrecher, Phys. Rev. B: Condens. Matter Mater. Phys., 2008, 77, 184417.
26. A.-P. Hynninen and M. Dijkstra, Phys. Rev. Lett., 2005, 94, 138303–138304.
27. A. Yethiraj, A. Wouterse, B. Groh and A. van Blaaderen, Phys. Rev. Lett., 2004, 92, 058301.
28. P. S. Colla, T. Mohanty, S. Nöjid, A. Riede, P. Schurtenberger and C. Likos, ACS Nano, 2018, 12, 4321–4337.
29. S. Semwal, C. Clowe-Coish, I. Saika-Voivod and A. Yethiraj, Phys. Rev. X, 2022, 12, 041021.
30. N. Li, H. Newman, M. Valera, I. Saika-Voivod and A. Yethiraj, Soft Matter, 2010, 6, 876–880.
31. E. A. Elfimova, A. O. Ivanov and P. J. Camp, J. Chem. Phys., 2012, 136, 194502.
32. H. R. Vutukuri, A. F. Demirörs, B. Peng, P. D. J. van Oostrum, A. Imhof and A. van Blaaderen, Angew. Chem., Int. Ed., 2012, 51, 11249–11253.
33. J. Wei, F. Sonmg and J. Dobiknar, Langmuir, 2016, 32, 9321–9328.
34. R. Messina and E. Kemgang, J. Chem. Phys., 2023, 159, 174903.
35. C. P. Royall and S. R. Williams, Phys. Rep., 2015, 560, 1–75.
36. I. A. Zaluzhnyy, R. P. Kurta, M. Scheele, F. Schreiber, B. I. Ostrovskii and I. A. Vartanyants, Mater. MDPI, 2019, 12, 3464.
37. A. van Blaaderen and P. Wiltzius, Science, 1995, 270, 1177–1179.
38. C. P. Royall, S. R. Williams, T. Ohtsuka and H. Tanaka, Nat. Mater., 2008, 7, 556.
39. M. Leocmach, J. Russo and H. Tanaka, J. Chem. Phys., 2013, 138, 12A515.
40. J. Taffs and C. P. Royall, Nat. Commun., 2016, 7, 13225.
41. F. Leoni and J. Russo, Phys. Rev. X, 2021, 11, 031006.
42. W. Gispen, G. M. Coli, R. van Damme, C. P. Royall and M. Dijkstra, ACS Nano, 2023, 17, 8807–8814.
43. J. Russo and H. Tanaka, Sci. Rep., 2012, 2, 505.
44. J.-P. Hansen and I. R. McDonald, Theory of Simple Liquids, Wiley, Chichester, 3rd edn, 2006.
45. J. F. Robinson, F. Turci, R. Roth and C. P. Royall, Phys. Rev. Lett., 2019, 122, 068004.
46. J. F. Robinson, F. Turci, R. Roth and C. P. Royall, Phys. Rev. E, 2019, 100, 062126.
47. C. Russ, K. Zahn and H. von Gruenberg, J. Phys.: Condens. Matter, 2003, 15, S3509.
48. C. P. Royall, J. Eggers, A. Furukawa and H. Tanaka, Phys. Rev. Lett., 2015, 114, 258302.
49. J. D. Honeycutt and H. C. Andersen, J. Phys. Chem., 1987, 91, 4950–4963.
50. M. Tanemura, Y. Hiwatari, H. Matsuda, T. Ogawa, N. Ogita and A. Ueda, Prog. Theor. Phys., 1977, 58, 1079–1095.
51. H. Jónsson and H. C. Andersen, Phys. Rev. Lett., 1988, 60, 2295–2298.
52. N. P. Bailey, J. Schiøtz and K. W. Jacobsen, Phys. Rev. B: Condens. Matter Mater. Phys., 2004, 69, 144205.
53. P. J. Steinhardt, D. R. Nelson and M. Ronchetti, Phys. Rev. B: Condens. Matter Mater. Phys., 1983, 28, 784–805.
54. S. Auer and D. Frenkel, J. Chem. Phys., 2004, 120, 3015–3029.
55. T. Kawasaki and H. Tanaka, Proc. Natl. Acad. Sci. U. S. A., 2010, 107, 14036–14041.
56. E. D. Cubuk, S. S. Schoenholz, J. M. Rieser, B. D. Malone, J. Rottler, D. J. Durian, E. Kaxiras and A. J. Liu, Phys. Rev. Lett., 2015, 114, 108001.
57. E. Boattini, S. Marin-Aguilar, S. Mitra, G. Foffi, F. Smallenburg and L. Filion, Nat. Commun., 2020, 11, 5479–5488.
58. E. Boattini, M. Dijkstra and L. Filion, J. Chem. Phys., 2019, 151, 154901.
59. G. Campos-Villalobos, E. Boattini, L. Filion and M. Dijkstra, J. Chem. Phys., 2021, 155, 174902.
60. J. Paret, R. L. Jack and D. Coslovich, J. Chem. Phys., 2020, 152, 144502.
61. F. C. Frank, Proc. R. Soc. A, 1952, 215, 43–46.
62. D. J. Wales, Energy Landscapes: Applications to Clusters, Biomolecules and Glasses, Cambridge University Press, Cambridge, 2004.
63. D. J. Wales and J. P. K. Doye, J. Phys. Chem. A, 1997, 101, 5111–5116.
64. L. Trombach, R. S. H. D. J. W. Hoy, L. Trombach, P. Schwerdtfeger, D. J. Wales and P. Schwerdtfeger, Phys. Rev. E, 2018, 97, 043309.
65. M. A. Miller and D. Wales, J. Phys. Chem. B, 2005, 109, 23109–23112.
66. M. A. Miller, J. J. Shepherd and D. J. Wales, Mol. Phys., 2008, 106, 1655–1664.
67. A. Malins, S. R. Williams, J. Eggers and C. P. Royall, J. Chem. Phys., 2013, 139, 234506.
68. K. Skipper, F. J. Moore and C. P. Royall, J. Chem. Phys., 2024, 2406, 16076.
69. J. Taffs, A. Malins, S. R. Williams and C. P. Royall, J. Phys.: Condens. Matter, 2010, 22, 104119.
70. C. L. Klix, K. Murata, H. Tanaka, S. R. Williams, A. Malins and C. P. Royall, Sci. Rep., 2013, 3, 2072.
71. J. Taffs, A. Malins, S. R. Williams and C. P. Royall, J. Chem. Phys., 2010, 133, 244901.
72. A. Malins, S. R. Williams, J. Eggers, H. Tanaka and C. P. Royall, J. Non-Cryst. Solids, 2010, 357, 760–766.
73. M. Godogna, A. Malins, S. R. Williams and C. P. Royall, Mol. Phys., 2010, 109, 1393–1402.
74. C. P. Royall, W. C. K. Poon and E. R. Weeks, Soft Matter, 2013, 9, 17–27.
75. C. P. Royall, S. R. Williams and H. Tanaka, J. Chem. Phys., 2018, 148, 044501.
76. A. I. Campbell and P. Bartlett, J. Colloid Interface Sci., 2002, 256, 325–330.
77. A. Hollingsworth, W. Russel, C. van Kats and A. van Blaaderen, APS March Meeting Abstracts, 2006, p. G21.005.
78. C. P. Royall, R. van Roij and A. van Blaaderen, J. Phys.: Condens. Matter, 2005, 17, 2315–2326.
79. A. Donovan, PhD thesis, University of Bristol, 2006.
80. W. C. K. Poon, E. R. Weeks and C. P. Royall, Soft Matter, 2012, 8, 21–30.
81. C. P. Royall, P. Charbonneau, M. Dijkstra, J. Russo, F. Smallenburg, T. Speck and C. Valeriani, Rev. Mod. Phys., 2024, 96, 045003.
82. J. Taffs, S. R. Williams, H. Tanaka and C. P. Royall, Soft Matter, 2012, 9, 297.
83. M. Leunissen, PhD thesis, Utrecht Universiteit, 2006.
84. A. Statt, R. Pinchaipat, F. Turci, R. Evans and C. P. Royall, J. Chem. Phys., 2016, 144, 144506.
85. J. C. Crocker and D. G. Grier, J. Colloid Interface Sci., 1995, 179, 298–310.
86. M. Leocmach and H. Tanaka, Soft Matter, 2013, 9, 1447–1457.
87. Y. Yang, nplocate, 2021, https://github.com/yangyushi/nplocate.
88. C. P. Royall, A. A. Louis and H. Tanaka, J. Chem. Phys., 2007, 127, 044507.
89. C. P. Royall and A. Malins, Faraday Discuss., 2012, 158, 301.
90. S. Plimpton, J. Comput. Phys., 1995, 117, 1–19.
91. A. Razali, C. J. Fullerton, F. Turci, J. E. Hallett, R. L. Jack and C. P. Royall, Soft Matter, 2017, 13, 3230.
92. R. N. Zia, B. J. Landrum and W. B. Russel, J. Rheol., 2014, 58, 1121–1157.
93. J. D. Weeks, D. Chandler and H. C. Andersen, J. Chem. Phys., 1971, 54, 5237.
94. J. A. Barker and D. Henderson, Rev. Mod. Phys., 1976, 48, 587–671.
95. J. P. K. Doye, D. J. Wales and R. S. Berry, J. Chem. Phys., 1995, 103, 4234–4249.
96. J. V. Selinger, Liquid Crystals, Springer International Publishing, Cham, 2016, pp. 131–182.
97. D. Coslovich, J. Chem. Phys., 2013, 138, 12A539.
98. A. Malins, S. R. Williams, J. Eggers, H. Tanaka and C. P. Royall, J. Phys.: Condens. Matter, 2009, 21, 425103.
99. T. M. Truskett, S. Torquato, S. Sastry, P. G. Debenedetti and F. H. Stillinger, Phys. Rev. E: Stat. Phys., Plasmas, Fluids, Relat. Interdiscip. Top., 1998, 58, 3083–3088.
100. C. P. Royall, A. Malins, A. J. Dunleavy and R. Pinney, J. Non-Cryst. Solids, 2015, 407, 34–43.
101. F. Leoni, F. Martelli and C. P. Royall, Phys. Rev. Lett., 2023, 130, 198201.
102. T. Jenkinson, P. Crowther, F. Turci and C. P. Royall, J. Chem. Phys., 2017, 147, 054501.
103. K. Thijssen, T. B. Liverpool, C. P. Royall and R. L. Jack, Soft Matter, 2023, 19, 7412–7428.
104. D. Coslovich, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2011, 83, 051505.
105. M. Leocmach and H. Tanaka, Nat. Commun., 2012, 3, 974.
106. A. Furukawa and H. Tanaka, Phys. Rev. Lett., 2010, 104, 245702.
107. J. de Graaf, W. C. K. Poon, M. J. Haughey and M. Hermes, Soft Matter, 2019, 15, 10–16.

---

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.

---