Coupling between criticality and gelation in "sticky" spheres: A structural analysis
David Richard James Hallett Thomas Speck C. Patrick Royall





Institut für Physik, Johannes Gutenberg-Universität Mainz, Staudingerweg 7-9, 55128 Mainz, Germany

HH Wills Physics Laboratory, Tyndall Avenue, Bristol BS8 1TL, UK

Institut für Physik, Johannes Gutenberg-Universität Mainz, Staudingerweg 7-9, 55128 Mainz, Germany

HH Wills Physics Laboratory, Tyndall Avenue, Bristol BS8 1TL, UK
School of Chemistry, University of Bristol, Cantock's Close, Bristol, UK
Centre for Nanoscience and Quantum Information, Tyndall Avenue, Bristol, UK


  We combine experiments and simulations to study the link between criticality and gelation in sticky spheres. We employ confocal microscopy to image colloid-polymer mixtures, and Monte Carlo simulations of the square-well (SW) potential as a reference model. To this end, we map our experimental samples onto the SW model. We find excellent structural agreement between experiments and simulations, both for locally favored structures at the single particle level and large-scale fluctuations at criticality. We follow in detail the rapid structural change of the critical fluid when approaching the gas-liquid binodal and highlight the role of critical density fluctuations for this structural crossover. Our results link the arrested spinodal decomposition to long-lived energetically favored structures, which grow even away from the binodal due to the critical scaling of the bulk correlation length and static susceptibility. 


Coupling between criticality and gelation in “sticky” spheres: A structural analysis
    C. Patrick Royall
    December 30, 2023
====================================================================================







§ INTRODUCTION


Understanding how an amorphous system becomes dynamically arrested upon compression or cooling is a long-standing challenge in statistical physics. Such amorphous solids encompass states of matter such as glasses, films, plastics, and gels, among others. Despite the fact that these systems are of technological importance and have received a lot of attention in the literature <cit.>, the microscopic mechanism responsible for macroscopic arrest remains elusive. For instance, the question whether the glass transition can be explained in the context of a thermodynamic or structural phase transition is still debated <cit.>. Interestingly, gels, in contrast to glasses, can have very sparse spatial structural arrangements, often described by a percolating network of bonded particles whose degree of dilution can in principle be unbounded <cit.>. These networks can result from cross-linking polymer chains or from physical bonds. The latter are caused, e.g., by depletion attractions <cit.> and for some time have served as a well controlled model system to study glasses and gels <cit.>.

When decreasing the range of attraction between particles, the gas-liquid coexistence becomes metastable with respect to fluid-solid coexistence, with a metastable critical point <cit.>. In the limit where the attraction range becomes smaller than roughly 10% of the diameter of the particles, the shape of the attractive part of the pair-wise potential becomes irrelevant <cit.>. This gave rise to the extended law of corresponding states <cit.>, which can also explain some aspects of the phase behavior of protein solutions <cit.>. Experimentally, such a short-ranged attraction can be probed using colloid-polymer mixtures <cit.>, so-called “sticky” spheres. Here, the polymers play the role of depletants for the large colloidal particles, where the radius of gyration of the polymer chains sets the range of attraction.

In the context of understanding the gelation of particles with short-ranged attraction, much work has been done during the last 20 years <cit.>. Several scenarios have been proposed, including diffusion-limited cluster aggregation <cit.>, phase separation <cit.>, and percolation <cit.>; and these can couple in different regions of the phase diagram <cit.>. Some of this work links the gelation to an arrested spinodal decomposition <cit.> and thus, a direct consequence of the underlying gas-liquid coexistence. A possible microscopic mechanism for the arrest was proposed with the presence, at gelation, of clusters that minimise the local potential energy <cit.>. Furthermore, it has been shown that gelation in sticky spheres offers a clear dynamical signature in comparison with a hard-sphere glass, where authors found a quasi-discontinuous increase in the relaxation time of the fluid for various packing fractions <cit.>. These dynamical transitions were found to be located at the gas-liquid binodal, even for a very dense gel, with a packing fraction exceeding the freezing point of the hard-sphere fluid <cit.>. Other work has suggested that the gelation line is located before the phase separation and extends at higher densities toward the location of the attractive glass <cit.>. A more recent numerical study <cit.> of the adhesive hard-sphere model <cit.> has linked the experimental gelation line to the mean-field rigidity transition introduced by He and Thorpe <cit.> in the context of random networks. More recently, a connection between directed percolation <cit.> as an equilibrium pre-structural transition to gelation has been proposed <cit.>.

The suggested relationship between the liquid-gas phase separation line <cit.> and gelation implies that critical fluctuations may influence the gelation process in the vicinity of the critical isochore. While criticality has been studied in colloid-polymer systems <cit.>, this work has tended to focus on systems where the interaction range is long enough that the system exhibits a stable colloidal liquid phase. In colloid-polymer mixtures, this corresponds to a polymer-colloid size ratio of around 0.3 <cit.>. For such systems, the colloidal liquid is not dense enough to arrest, and gelation is only achieved upon quenching with a much stronger attraction strength than that requited for criticality <cit.>. However, in the sticky sphere limit, immediately upon quenching through the binodal, the density of the colloidal “liquid” is sufficient that the system undergoes dynamical arrest <cit.>. Under these circumstances, gelation may couple to critical fluctuations and this forms the subject of our study. 





§ METHODS





 §.§ Experiments


We employ confocal microscopy and particle tracking to resolve the positions of the colloidal particles. Colloid-polymer mixtures are composed of polystyrene (PS) polymer chains and sterically stabilized polymethyl methacrylate (PMMA) spheres with diameter σ=2950 nm and polydispersity Δ=5%, determined via scanning electronic microscopy (SEM). We use rhodamine as fluorescent label. Polystyrene has a molecular weight M_w=1.3×10^6 corresponding to an effective radius of gyration R_g=35 nm under θ conditions. We use a solvent mixture of cis-decalin and cyclohexyl bromide which is density and refractive index-matched. We additionally screen electrostatic interactions using 4 mMol of tetrabutyl ammonium bromide salt. From previous work <cit.>, we estimate R_g≃50 nm at room temperature. This leads to a polymer-colloid size ratio q≃0.03 approaching the sticky sphere limit. We then seal each sample into a borosilicate glass capillary with epoxy resin.

Samples are imaged by confocal microscopy using a Leica SP8. We image different parts of the suspension at least 15 μm from the wall. To extract the colloidal particle positions we employ three dimensional particle tracking using the difference of Gaussian method. This tracking is performed using the package ”Colloid” developed by Leocmach et al. <cit.>. From the tracking we estimate the sample packing fraction as ϕ≃π Nσ^3/(6V), with N the number of tracked particles and V the volume of the sample. It is well known that sterically stabilized PMMA colloidal particles can exhibit swelling and unscreened electrostatic interactions <cit.> which can lead to an intrinsic softness modeled by an effective diameter σ_eff>σ. To estimate σ_eff under our experimental conditions, we apply a method similar to the one used in Ref. <cit.>. In this study, the authors matched the pair correlation function g(r) of dense hard spheres with the known Percus-Yevick expression by varying the effective diameter σ_eff. In our work, we extend this mapping with matching the g(r) of the colloid-polymer mixture with an attractive square-well fluid, which we adopt throughout this study as a reference system for our mixtures. In addition, our mapping allows us to determine the systematic tracking errors of the colloidal positions, responsible for the broadness of g(r) at contact <cit.>. More details can be found in the appendix <ref>. We find σ_eff=3100  nm, which corresponds to an effective diameter that is ∼5% larger than for the SEM estimation.



 §.§ Computer simulations


We study the behavior of the square-well (SW) model serving as a reference model for our colloid-polymer mixtures. We perform standard Monte Carlo simulations in the NVT ensemble employing local moves. The system is composed of N=5000 particles with diameters drawn from a Gaussian distribution with polydispersity Δ=5%. The interaction potential between two particles i and j is


    V(r)=∞   if r ≤σ_ij
    
    -V_SW   if σ_ij < r <σ_ij+δ
    
    0    if r ≥σ_ij+δ,


where σ_ij=(σ_i+σ_j)/2, σ_i and σ_j being the diameter of particle i and j respectively. The attraction range δ is set by the polymers and fixed to 0.03σ and independ of σ_ij. If not mentioned otherwise, simulations are performed in a cubic box. Additionally, we use a slab geometry to compute the coexistence packing fraction between the gas and liquid phase with box lengths: L_x=L_y=L_z/2. We determine the binodal by fitting the density profile along the gas-liquid coexistence through,


    ϕ(z)=ϕ_l-ϕ_g/2+ϕ_l-ϕ_g/2tanh(z-z_0/2w).


Here, ϕ_g and ϕ_l denote the gas and liquid coexistence packing fractions, and z_0 and w the interface position and width, respectively. Since for this polydispersity the liquid is still metastable with respect to the crystal, we check that all liquid slabs remained in the liquid phase. All distances and energies are expressed in units of σ and k_BT, k_B being the Boltzmann constant. Throughout, we denote the dimensionless attraction strength by V=V_SW/k_BT.



 §.§ Mapping procedure


The attraction strength between colloidal particles is controlled by the polymer concentration c_p. However, it is challenging to determine c_p precisely enough to map it directly. Instead, we map every sample individually through matching the experimentally measured total correlation function h_exp(r)=g_exp(r)-1 and distribution P_exp(n) of bond number n to the SW fluid <cit.>, yielding an effective attraction strength V for every experimental sample.



Specifically, we compute rh_sim(r) and  P_sim(n) on a grid in the (ϕ, V) plane around the critical point. We evaluate P(n) by constructing a bond network using a Voronoi decomposition considering only direct Voronoi neighbors <cit.> with the bond distance cutoff r_c set to 1.5σ. Following the work of Largo et al., we estimate the location of the critical point to be at ϕ_c≃0.275 and V_c≃3.22 for zero polydispersity <cit.>. We shall show later on that this value is very close to the critical point of our model with Δ=5%. We use packing fractions ranging from ϕ=0.05 to 0.5 with an interval Δϕ=0.025 and V values ranging from 0 to 4 with Δ V=0.2 for 0<V<2 and Δ V=0.1 for 2<V<4. Overall, we end up with a grid of 608 state points. Beyond V≃3.2, the fluid crosses the binodal and starts to form a gel. Thus, structural observables such as h(r) and P(n) might evolve due to aging. We fix for every state point a MC relaxation time of 10^6 steps before we compute any observables. Another 10^6 steps is used to compute observables. We then pick for each sample the two numerical packing fractions  ϕ_- and ϕ_+, which encompass our sample density. We then compute for each V value h_sim and P_sim as a linear combination of ϕ_- and ϕ_+. We finally compute as a “goodness” parameter for our matching procedure the least squares


    χ^2_h = ∑_i [ r_i(h_exp(r_i) - h_sim(r_i))]^2

and

    χ^2_P = ∑_i [P_exp(n_i) - P_sim(n_i)]^2.


In practice, χ^2_h is computed for r<4σ, while an additional Gaussian (Δ=5%) noise is added on numerical positions to mimic particle tracking errors <cit.>, see appendix <ref>. The global minima of χ^2_h and χ^2_P give us two independent evaluations of V_exp. We then assign for each sample the mean of those two, V_exp=[V_h+V_P]/2, whereas (|V_h-V_P|)/2 serves as an estimation for errors.

In Fig. <ref>(a) and (b), we present χ^2_h and χ^2_P as functions of V for three samples picked along the critical isochore. From a hard sphere fluid (gray curves) to the gel (red curves), we always observe a global minimum for χ^2 and find a good agreement between V_h and V_P. Results of the SW model are shown in Fig <ref>(c,d) and compared to experimental data. The correspondence between experiments and simulations is excellent, we observe only small discrepancies in the bond-number distribution for the gel sample for small bond numbers n, see Fig <ref>(d). These deviations are not surprising since we do not include the effect of aging in our mapping method.





§ RESULTS





 §.§ Phase diagram





We start by discussing the phase behavior of our system in Fig. <ref>. We choose to present our data in the ϕ-V plane. The gas-liquid binodal of our reference system is indicated by the black solid line. The location of the critical point is indicated by a black star and was determined via the block distribution functions method <cit.>, we shall come back to this point later. From our mapping procedure, we can place each sample at a given ϕ and V. We distinguish an ergodic fluid in contrast to a gel by respectively blue and red square symbols. We also indicate samples without polymers, i.e. hard spheres, in gray. The typical gel sample shows a clear dynamical arrest with a relaxation time beyond 100 τ_B (Supplementary Video) with τ_B  being the Brownian time. This is consistent with a previous dynamical experimental study for a very similar system <cit.>.

The main observation here is that arrested samples are lying close to the binodal for a wide range of densities (0.1<ϕ<0.4). That is to say, all samples which are gels always map onto state points of the reference system which are on or above the binodal. This means that the radial distribution function and bond distribution of a gel sample never correspond to an ergodic fluid. This is in agreement with previous numerical and experimental studies, which associate the gelation to the location of the spinodal and thus a direct consequence of equilibrium properties of the gas-liquid coexistence <cit.>. For short-ranged attractive systems, the binodal and spinodal are located very close together, such that in experiments it is hard to distinguish them <cit.>. The highest attraction strength at which the sample remains ergodic is found at V≃2.9. Additionally, in Fig. <ref>, we show images from confocal microscopy along the critical isochore ϕ≃ϕ_c. We observe a rather continuous change of structural behavior from hard-sphere fluid to gel. At V=0 [Fig. <ref>(b)] the spatial distribution of colloids is homogeneous at length scales larger than the particle size. When increasing the concentration of polymers, i.e. increasing V, clusters start to form with lower density in between [Fig. <ref>(c-d)]. Finally, the suspension shows an arrested spinodal decomposition, which is revealed in Fig. <ref>(e) by large assemblies rich in colloidal particles and depleted zones without any particles.



 §.§ Critical density fluctuations



In section <ref>, we have shown that the gelation line is located close to the binodal and that dense domains of colloidal particles grow progressively as a function of the polymer concentration. We now demonstrate how this change of behavior can be directly linked to the fluctuations present in the context of criticality. To quantify the spatial evolution of the density we use the block distribution functions method <cit.>, which provides the location of the critical point V_c <cit.> and to some extent the isothermal susceptibility χ of the fluid <cit.>. The procedure is as follows: We divide our system into a series of cubic subcells of dimension l=L/b, with b being an integer. The global density is defined as ρ = 1/n_b∑_i ρ_i, where n_b=b^3 is the number of subcells and ρ_i is the local density in subcell i. We then extract the second and forth moment of the density ⟨ m_2 ⟩ and ⟨ m_4 ⟩ computed as ⟨ m_x ⟩_l = 1/n_b∑_i (ρ_i-ρ)^x. We can finally define the Binder cumulant through


    U_l = ⟨ m_4 ⟩_l/⟨ m_2 ⟩_l^2.


One of the main properties of U_l is its size invariance with respect to l at the critical point <cit.>. This allows us to determine unambiguously and accurately the critical attraction strength V_c in computer simulations. In more detail, we compute U_l for various attraction strengths V and subcell lengths l. To this end, we sample U_l for each V values at a fixed packing fraction ϕ≃ϕ_c=0.275 using an extra 10^8 MC steps. In Fig. <ref>, we show the evolution of U_l as a function of V. We observe a crossing for V>3, which can be resolved more accurately as shown in the inset of Fig. <ref>(a). We find V_c=3.19(3) as the final value for the critical point, which is close to the value determined by Largo et al. <cit.> for zero polydispersity. Additionally, in Fig. <ref> we show experimental results for U_l at a fixed subcell length, l=2σ. Although the experimental data suffer from a lack of statistic in comparison with simulation data, we observe an overall good agreement. More precisely, gel samples indicated by red squares are, within the errors, located either at the cumulant crossing or at larger attraction strengths.



In practice, the change of U_l indicates that the density distribution P_l(ρ_i) moves from a Gaussian shape, centered at ρ_i=ρ, to a binodal shape, where the two maxima of the distribution move progressively towards the coexistence densities ρ_- and ρ_+. We show such evolution for both experiments and simulations in Fig. <ref>(a-d) for subcells of length l≃3σ.  To be consistent with the phase diagram in Fig. <ref>, we choose to plot P_l(ϕ) instead of P_l(ρ). We start with a hard-sphere sample [Fig. <ref>(a)] for which we observe a narrow distribution lacking low (ϕ<0.1) and high density (ϕ>0.5) regions. The same behavior continues up to V=2 [Fig. <ref>(b)]. When further increasing V, we start to observe broader distributions with almost empty (ϕ<0.1) and colloid-rich regions (ϕ>0.5). This can be clearly seen for our last ergodic sample at V=2.9 [Fig. <ref>(c)]. Interestingly, dense regions can reach packing fractions larger than the freezing point of a hard-sphere fluid. This explains why for non-polydisperse samples, one can expect a speed up of crystallization around the critical point, where the dense regions arising from criticality will lower the nucleation barrier <cit.>. Continuing to quench, the distribution becomes broader and non-Gaussian at V=3.2, which also corresponds to the location of our first gel sample. Surprisingly, we find that P_l^exp matches quite closely the equilibrated P_l^sim, indicating that aging effects are not pronounced. Finally, we can say that the overall behavior of P_l(ϕ) confirms quantitatively the observation made for confocal images in the previous section [Fig. <ref>(a-d)].





 §.§ Fractal dimension, bulk correlation length, and bond distribution




We now turn to discuss the characteristic cluster shape and length scale formed by colloidal particles when quenched though criticality. We first compute the scattering intensity I(q) as a function of the wave vector q from the Fourier transform of the confocal pixel map. In Fig. <ref>(a) we show a typical colloidal gel close to criticality and its associated Fourier spectrum in the inset. A radial average leads to I(q) as show in Fig (b) for different samples along the critical isochore. When quenching, we observe a divergence of I(q) for small q as a consequence of larger domains forming and diffusing light. Since the scattering intensity I(q) is proportional to P(q)S(q), where P(q) is the form factor and S(q) the structure factor. We expect for q>2π/σ the scaling I(q)∼ q^-4, which is confirmed by the dashed line in Fig. <ref>(b). For small q, the scaling depends on the shape of the diffusing clusters through d_f, its the fractal dimension. Previous small-angle neutron scattering (SANS) experiments have shown that the fractal dimension at the gelation is typically around d_f≃1.7 <cit.>, which is supported by particle-resolved experiments <cit.> and simulation <cit.>.

We will discuss this point in detail later on. We can additionally extract S(q) for small wave vector from S(q)=I(q)/(AP(q)), where we evaluate AP(q)) by fitting I(q) for q>2π/σ by a polydisperse form factor and a prefactor A, as shown in the inset of Fig. <ref>(b). The resulting procedure is shown in Fig (c) and directly compared to simulation data. We find an overall very good agreement between experimental data and simulations. The structure factor can then be used to extract an estimation of the bulk correlation length ξ through the Orstein-Zernike scaling,

    S(q)=S_0/1+(ξ q)^2,

which holds close to criticality for small q. In Fig. <ref>(d), we plot the evolution of the fractal dimension d_f as a function of the attraction strength V. We find when increasing V a progressive increase of d_f and saturates for gel samples at d_f=2.4-2.5. These values agree with fractal dimensions found for gel at low volume fractions and large attraction strength <cit.>. At criticality, where gelation occurs, we find d_f=1.6-1.7, which is consistent with a previous SANS experimental study <cit.>.

In Fig. <ref>(e), we present the behavior of the bulk correlation length ξ as a function of the reduced attraction strength ϵ. We find for simulation data that ξ diverges for ξ→0 and that values are well modelled by the Ising universality class where ξ∼ϵ^-ν up to ϵ≃0.2, where in 3D ν≃0.63. In contrast, for the experimental data the equilibrium fluids follow closely the simulation data until gelation occurs close to the binodal, whereby the system is arrested and therefore the correlation length does not increase beyond ξ≈2σ. This rather small spatial correlation is not in contradiction with confocal images of gels, where structures are quite ramified in a network without large colloidal domains, cf. Fig. <ref>(e) and Fig. <ref>(a). This behavior would indicates that gels are, in a way, pictures of early critical fluids, where the cost of breaking bonds does prevent relaxation and thus the growth of the correlation length ξ.

Finally, in Fig. <ref>(d) we discuss the behavior of the bond distribution through its variance ⟨ n^2⟩. It was shown recently for a similar colloid-polymer mixture model based on the Asakura-Oosawa potential that the variance of the bond distribution as a function of V is peaked close to gelation <cit.>. We observe the same behavior in our simulations crossing the binodal, where ⟨ n^2⟩ exhibits a maximum for V≃3.5. Quenching further the variance decreases, indicating a progressive aging of the network structure. The behavior of the variance will of course depend on time and other factors such as hydrodynamic interactions <cit.>. In experiments, we also observe a growth of ⟨ n^2⟩ for equilibrium fluids (V<3). We find larger variances for gel samples, but not as high as in simulations, which also confirms the picture found for ξ in Fig. <ref>(e), where the gel structure corresponds to an arrested critical fluid. We return to the behaviour in the gel state in Sec. <ref>.



 §.§ Local structure





We now go to smaller scales and inquire what type of local structures emerge in a fluid close to criticality, and on their role in the dynamical arrest of the gel network. We have seen that in our experiments the typical length scale before gelation reaches ≈2σ, emphasizing that the change in the dynamics at gelation is quite local and might only be caused by the first neighbor shell of each particle. To gain insight into the local structure of the super-critical fluid, we employ the topological cluster classification (TCC) <cit.>. We use the same bond network as described in the methods section. From the TCC, we extract the population of clusters of size m composed of 3 to 13 particles. This includes crystal structures like face-centered cubic (fcc) and hexagonal close-packed (hcp). In the TCC, each cluster type is labeled with two characters mX, with m the number of particles composing the cluster and X a letter to distinguish different spatial symmetries. Those letters follow the convention of the minimum energies of various potentials such as the Morse potential <cit.>.  We then identify for each cluster size m the predominant type X which corresponds to a minimum energy structure, a locally favored structure (LFS) <cit.>. In practice, we determine these clusters by following the global population spectrum of all clusters at different values of the attraction strength V. The population frequency of a cluster of type c is computed as N_c/N, with N_c the number of particles composing a cluster of type c and N the total number of particles. In the experimental data, we do not consider colloidal particles that are less than σ away from an edge.

In Fig. <ref>(a), we show the population of clusters comparing experiments and simulations. We find an excellent agreement, which demonstrates that mapping to the SW fluid reproduces also the higher-order structural features present in the experiments. We pick three different state points along the critical isochore including a hard-sphere fluid (V=0), near-critical fluid (V=2.9), and gel (V=3.2). We first notice that the overall histograms show that when increasing V (from top to bottom), larger clusters composed of more than 10 particles appear. Additionally, we do not observe any signature of crystallization, i.e. no fcc nor hcp. However, we find that medium size clusters with 5<m<10 grow progressively with respect to V, and always with the same order among clusters with equal number m of particles. In the following, we exclude larger clusters with m>12 from the discussion since their weights are negligible.



§ DISCUSSION





 §.§ Entropy favors low-symmetry clusters


For isolated clusters it has been shown that the relative population of clusters with the same number of bonds is determined by entropy <cit.>. Interestingly, this still holds for the fluid and even the gel as shown in Fig. <ref>(b) for m=6. We find that 6Z clusters, which are polytetrahedra, are always significantly more prevalent in the system compared to 6A clusters, which are octahedrons (the population goes up from ≃5% for isolated 6A to ≃10% in the fluid). The same trend can be observed for other clusters, with the more symmetric clusters being less populated. This demonstrates that the minimum (free) energy clusters are determined by rotational and vibrational entropies. It has also been shown that the potential part of the free energy will promote both octahedral and tetrahedral order. The same observation was found in a more recent numerical study of the SW fluid <cit.>, where gelation was associated with polytetrahedral order. Another important observation is that we also find a large population of clusters with fivefold symmetry (8B, 9B, and 10B) at any V. These structures are known to be local energy minima of the Morse potential <cit.>. Hence, even though a gel is in an energy landscape far away from the equilibrium state, it can still locally minimize its free energy by forming isolated locally favored structures, which eventually will overlap and form the gel network. For larger clusters with 11 and 12 particles, we observe several Morse minima that correspond to 10<ρ_0<25 (ρ_0 is the Morse potential parameter controlling the range of attraction, in our case δ). The range for ρ_0 found here is consistent with short range attraction. We found that 11E and 12D are the dominant structures for m=11 and m=12 respectively, which corresponds to ρ_0≃17 <cit.>.



 §.§ Hard sphere fluid


Having obtained a set of minimum energy clusters, we use these to follow the overall structural change in the fluid. As a reference, we first consider the compression of the hard sphere fluid (V=0) towards the glass shown in Fig. <ref>(a). At low density, ϕ≃0.05, there is an absence of clusters as one can expect since particles are mostly isolated without any neighbors. At this low density, the system exhibits as a larger cluster only a few percents of 4A (tetrahedron). The typical spatial arrangement of clusters is illustrated with an experimental snapshot in Fig. <ref>(a1). We then observe a rapid change from ϕ=0.05 to 0.25, where medium sized clusters such as 5A, 6Z, 7C, and 8B start to appear, cf. Fig. <ref>(a1-a3). This continues until a crossover at ϕ>0.25, where isolated and small clusters (3A and 4A) are converted with a combination of medium clusters (5A, 6Z, 7C, and 8B) to larger structures: 9B, 10B, 11E, and 12D, see Fig. <ref>(a4). This kind of conversion will continue at higher packing fractions with the extinction of 5A, 6Z, and 7C to promote clusters sharing the same sub-structures. We also notice that the agreement between experiments and simulation is excellent. This gives another confirmation for the robustness of the mapping procedure employed.







 §.§ En route to the gel



In Fig. <ref>(b), we follow the same idea but now approaching the gel. We fix the density to the critical packing fraction ϕ_c=0.275 and progressively increase the attraction strength V towards the binodal. We find a plateau with little change in structure until V≃2. We show the spatial arrangement of clusters at V=1.1 in Fig. <ref>(b1), where we only display clusters with m>4. For higher attraction strengths, V>2, we observe a structural crossover (cf. the compression of hard spheres). Clusters with m>4 increase quickly until reaching the binodal. Interestingly, this change starts where we have located the beginning of the critical scaling for the bulk correlation length ξ and the static susceptibility χ. The main contribution to this structural evolution comes from the creation of a large amount of 6Z, 7C, and 8B clusters, see Fig. <ref>(b1-b3). These clusters are not necessarily localized in space, but we do observe denser regions with even the presence of 9B at V=2.9, see Fig. <ref>(b3). This is consistent with our previous finding in Fig. <ref>(c), where the density distribution exhibits a broader tail toward large packing fractions ϕ>0.4. Therefore, we can expect a link between the local structures of dense hard spheres and the colloidal-rich regions arising from critical density fluctuations. Finally at V>V_c, which corresponds to our first identification of a gel, clusters percolate the whole system, see Fig. <ref>(b4). The inner parts of the percolating structure are rich in large clusters such as 8B, 9B, and 10B, which are defective icosahedra sharing a fivefold symmetry. They are known to play an important role in the slowing down of the dynamics for dense hard-spheres approaching the glass transition <cit.>. This may help to explain the origin of the rigidity of the network leading to gelation <cit.>.

The crucial difference between a gel and a dense hard sphere liquid is that not only the defective icosahedra (8B, 9B, and 10B) are very stable due to entropy, but each bond breaking will also result in an energetic penalty of more than 3k_BT. Therefore only colloids in the outer part of dense regions will be able to break bonds and diffuse (see Supplementary Video). To conclude, we find that the structural crossover starts 1k_BT before the location of the binodal when increasing V. This is a direct consequence of criticality inducing larger density fluctuations. For a monodisperse sample, these fluctuations might promote crystal precursors and lower the nucleation time <cit.>. For a polydisperse sample, however, they induce low-symmetry polytetrahedral backbones for the bond network. These clusters are incompatible with respect to the crystal symmetry. Hence, the system will fall into an amorphous solid state, the gel. Regarding the mapping we comment that, since the fluid is ergodic, the mapping to the SW fluid should be robust (as indeed it is). In contrast, the gel is non-ergodic and thus the correspondence would depend on the dynamics and history <cit.>.



§ CONCLUSIONS



We have investigated the role of criticality in the gelation of sticky spheres. Combining experiment and simulations, we provide further evidence that the dynamical arrest is initiated by the onset of critical fluctuations in agreement with previous work <cit.>. We have demonstrated that carefully mapping two-point structure to the square-well fluid faithfully reproduces the experimental data including higher-order local structures as identified by the topological cluster classification method <cit.>. In particular, we find that gel samples are (i) located at the cumulant crossing, (ii) identified by a broad distribution of densities, and (iii) have correlation lengths of ξ≃2σ. We find a sharp but continuous increase of locally favored structures when increasing the attraction strength. This increase occurs in concert with the increase of both the bulk correlation length and the static susceptibility, which can be extracted together with the structure factor of the fluid. More precisely, the start of the critical scaling of these two quantities coincides with the appearence of larger locally favored structures. The picture of a gel is thus that of an early critical fluid, which is arrested due to the large cost of breaking bonds. Before arrest, clusters of several particles appear, which have a low symmetry favored by entropy. The densification of these clusters driven by the incipient critical fluctuations then leads to the gel.





We gratefully acknowledge Francesco Turci and Matteo Campo for stimulating discussions. We thank M. Schmiedeberg for helpful discussions. We acknowledge financial support by the DFG through collaborative research center TRR 146 (D.R.).





§ EFFECTIVE COLLOID DIAMETER



There are two challenges for the determination of the effective packing fraction and the effective temperature of a colloidal sample from real space imaging. First, it is difficult to know the change of the colloidal diameter in solution due to swelling and unscreened electrostatic <cit.>. This leads to a poor estimation of the sample packing fraction and, e.g., a mismatch of structural oscillations in pair correlations. Second, errors due to imaging and tracking will lead to an error on the true position of a particle, resulting in a broadening of the peaks of a pair correlation <cit.>. In the context of matching a sample's g(r) onto simulations, these errors will induce a systematic overestimation of the temperature of a mixture at a given polymer concentration c_p. To handle these issues, we pick a sample close to criticality and compute its pair correlation r[g(r)-1] for various effective diameters. To mimic tracking errors, we pick a simulation state point at criticality and apply different Gaussian noises of variance Δ_error to the particle positions. We then find the optimal set of parameters that minimize χ^2 = ∑_i [ r_i(h_exp(r_i) - h_sim(r_i))]^2. The evolution of χ^2 as a function of the parameters is shown in Fig. <ref>(a). We find a unique minimum for χ^2 leading to an effective diameter σ=3.1μ m (5% larger than for dry colloids) and a tracking error Δ_error=5%, which is consistent with previous works <cit.>. The resulting matched pair correlation function is shown in Fig. <ref>(b). We observe a very good agreement between experiment and simulation which validate σ and Δ_error.





69
natexlab#1#1bibnamefont#1#1bibfnamefont#1#1citenamefont#1#1url<#>1urlprefixURL 


[Debenedetti and
  Stillinger(2001)]debenedetti2001supercooled
authorP. G. Debenedetti
  and authorF. H.
  Stillinger, journalNature
  volume410, pages259 (year2001).

[Berthier and Biroli(2011)]berthier2011theoretical
authorL. Berthier and
  authorG. Biroli,
  journalRev. Mod. Phys. volume83,
  pages587 (year2011).

[Biroli and Garrahan(2013)]biroli2013perspective
authorG. Biroli and
  authorJ. P. Garrahan,
  journalJ. Chem. Phys volume138,
  pages12A301 (year2013).

[Berthier and Ediger(2016)]berthier2016
authorL. Berthier and
  authorM. D. Ediger,
  journalPhys. Today volume69,
  pages40 (year2016).

[Turci et al.(2017)Turci, Royall, and
  Speck]turci2017prx
authorF. Turci,
  authorC. P. Royall,
  and authorT. Speck,
  journalPhys. Rev. X volume7,
  pages031028 (year2017).

[Royall and Williams(2015)]royall2015role
authorC. P. Royall and
  authorS. R. Williams,
  journalPhys. Rep volume560,
  pages1 (year2015).

[Royall et al.(2017)Royall, Turci,
  Tatsumi, Russo, and Robinson]royall2017
authorC. P. Royall,
  authorF. Turci,
  authorS. Tatsumi,
  authorJ. Russo, and
  authorJ. F. E. Robinson,
  journalArXiV p. pages1711.04739
  (year2017).

[Manley et al.(2005)Manley, Wyss,
  Miyazaki, Conrad, Trappe, Kaufman, Reichman, and
  Weitz]manley2005glasslike
authorS. Manley,
  authorH. Wyss,
  authorK. Miyazaki,
  authorJ. Conrad,
  authorV. Trappe,
  authorL. Kaufman,
  authorD. Reichman, and
  authorD. Weitz,
  journalPhys. Rev. Lett. volume95,
  pages238302 (year2005).

[Griffiths et al.(2017)Griffiths, Turci,
  and Royall]griffiths2017local
authorS. Griffiths,
  authorF. Turci, and
  authorC. P. Royall,
  journalJ. Chem. Phys volume146,
  pages014905 (year2017).

[Poon(2002a)]poon2002physics
authorW. Poon, journalJ.
  Phys. Condens. Matter volume14, pagesR859
  (year2002a).

[Van Megen and Pusey(1991)]van1991dynamic
authorW. Van Megen and
  authorP. Pusey,
  journal‎Phys. Rev. A volume43,
  pages5429 (year1991).

[De Hoog et al.(2001)De Hoog, Kegel,
  Van Blaaderen, and Lekkerkerker]de2001direct
authorE. De Hoog,
  authorW. Kegel,
  authorA. Van Blaaderen,
  and
  authorH. Lekkerkerker,
  journalPhys. Rev. E volume64,
  pages021407 (year2001).

[Cipelletti and Ramos(2005)]cipeletti2005
authorL. Cipelletti and
  authorL. Ramos, journalJ.
  Phys.: Condens. Matter volume17,
  pagesR253–R285 (year2005).

[Zaccarelli(2007a)]zaccarelli2007colloidal
authorE. Zaccarelli,
  journalJ. Phys. Condens. Matter volume19,
  pages323101 (year2007a).

[Meng et al.(2010)Meng, Arkus, Brenner,
  and Manoharan]meng2010free
authorG. Meng,
  authorN. Arkus,
  authorM. P. Brenner,
  and authorV. N.
  Manoharan, journalScience
  volume327, pages560 (year2010).

[Royall et al.(2008)Royall, Williams,
  Ohtsuka, and Tanaka]royall2008direct
authorC. P. Royall,
  authorS. R. Williams,
  authorT. Ohtsuka, and
  authorH. Tanaka,
  journalNat. Mater volume7,
  pages556 (year2008).

[Taylor et al.(2012)Taylor, Evans, and
  Royall]taylor2012temperature
authorS. L. Taylor,
  authorR. Evans, and
  authorC. P. Royall,
  journalJ. Phys. Condens. Matter volume24,
  pages464128 (year2012).

[Pinchaipat et al.(2017)Pinchaipat,
  Campo, Turci, Hallet, Speck, and Royall]pinchaipat2017
authorR. Pinchaipat,
  authorM. Campo,
  authorF. Turci,
  authorJ. Hallet,
  authorT. Speck, and
  authorC. P. Royall,
  journalPhys. Rev. Lett. volume119,
  pages028004 (year2017).

[Haxton et al.(2015)Haxton, Hedges, and
  Whitelam]haxton2015crystallization
authorT. K. Haxton,
  authorL. O. Hedges,
  and authorS. Whitelam,
  journalSoft matter volume11,
  pages9307 (year2015).

[Largo et al.(2008)Largo, Miller, and
  Sciortino]largo2008vanishing
authorJ. Largo,
  authorM. Miller, and
  authorF. Sciortino,
  journalJ. Chem. Phys volume128,
  pages134513 (year2008).

[Noro and Frenkel(2000)]noro2000extended
authorM. G. Noro and
  authorD. Frenkel,
  journalJ. Chem. Phys volume113,
  pages2941 (year2000).

[Stradner et al.(2004)Stradner, Sedgwick,
  Cardinaux, Poon, Egelhaaf, and Schurtenberger]stradner2004equilibrium
authorA. Stradner,
  authorH. Sedgwick,
  authorF. Cardinaux,
  authorW. C. Poon,
  authorS. U. Egelhaaf,
  and
  authorP. Schurtenberger,
  journalNature volume432,
  pages492 (year2004).

[Katsonis et al.(2006)Katsonis, Brandon,
  and Vekilov]katsonis2006corresponding
authorP. Katsonis,
  authorS. Brandon, and
  authorP. G. Vekilov,
  journalJ. Phys. Chem. B volume110,
  pages17638 (year2006).

[Platten et al.(2015)Platten,
  Valadez-Pérez, Castañeda-Priego, and Egelhaaf]platten2015extended
authorF. Platten,
  authorN. E. Valadez-Pérez,
  authorR. Castañeda-Priego,
  and authorS. U.
  Egelhaaf, journalJ. Chem. Phys
  volume142, pages05B602_1
  (year2015).

[Grant and Russel(1993)]grant1993volume
authorM. Grant and
  authorW. Russel,
  journalPhys. Rev. E volume47,
  pages2606 (year1993).

[Verduin and Dhont(1995)]verduin1995phase
authorH. Verduin and
  authorJ. K. Dhont,
  journalJ. Colloid Interface Sci volume172,
  pages425 (year1995).

[Poon(2002b)]poon2002
authorW. C. K. Poon,
  journalJ. Phys.: Condens. Matter volume14,
  pagesR859 (year2002b).

[Zaccarelli(2007b)]zaccarelli2007
authorE. Zaccarelli,
  journalJ. Phys.: Condens. Matter volume19,
  pages323101 (year2007b).

[Lu et al.(2008)Lu, Zaccarelli, Ciulla,
  Schofield, Sciortino, and Weitz]lu2008gelation
authorP. J. Lu,
  authorE. Zaccarelli,
  authorF. Ciulla,
  authorA. B. Schofield,
  authorF. Sciortino,
  and authorD. A. Weitz,
  journalNature volume453,
  pages499 (year2008).

[Laurati et al.(2009)Laurati, Petekidis,
  Koumakis, Cardinaux, Schofield, Brader, Fuchs, and Egelhaaf]laurati2009
authorM. Laurati,
  authorG. Petekidis,
  authorN. Koumakis,
  authorF. Cardinaux,
  authorA. B. Schofield,
  authorJ. M. Brader,
  authorM. Fuchs, and
  authorS. U. Egelhaaf,
  journalJ. Chem. Phys. volume130,
  pages134907 (year2009).

[Eberle et al.(2011)Eberle, Wagner, and
  Castañeda-Priego]eberle2011dynamical
authorA. P. Eberle,
  authorN. J. Wagner,
  and
  authorR. Castañeda-Priego,
  journalPhys. Rev. Lett. volume106,
  pages105704 (year2011).

[Pandey and Conrad(2013)]pandey2013
authorR. Pandey and
  authorJ. C. Conrad,
  journalSoft Matter volume9,
  pages10617 (year2013).

[Royall et al.(2015a)Royall,
  Eggers, Furukawa, and Tanaka]royall2015probing
authorC. P. Royall,
  authorJ. Eggers,
  authorA. Furukawa, and
  authorH. Tanaka,
  journalPhys. Rev. Lett. volume114,
  pages258302 (year2015a).

[Royall et al.(2018)Royall, Williams, and
  Tanaka]royall2018vitri
authorC. P. Royall,
  authorS. R. Williams,
  and authorH. Tanaka,
  journalJ. Chem. Phys volume148,
  pages044501 (year2018).

[Bibette et al.(1992)Bibette, Mason,
  Gang, and Weitz]bibette1992kinetically
authorJ. Bibette,
  authorT. Mason,
  authorH. Gang, and
  authorD. Weitz,
  journalPhys. Rev. Lett. volume69,
  pages981 (year1992).

[Kroy et al.(2004)Kroy, Cates, and
  Poon]kroy2004cluster
authorK. Kroy,
  authorM. Cates, and
  authorW. Poon,
  journalPhys. Rev. Lett. volume92,
  pages148302 (year2004).

[Valadez-Pérez
  et al.(2013)Valadez-Pérez, Liu, Eberle, Wagner, and
  Castañeda-Priego]valadez2013dynamical
authorN. E. Valadez-Pérez,
  authorY. Liu,
  authorA. P. Eberle,
  authorN. J. Wagner,
  and
  authorR. Castañeda-Priego,
  journalPhys. Rev. E volume88,
  pages060302 (year2013).

[Kohl et al.(2016)Kohl, Capellmann,
  Laurati, Egelhaaf, and Schmiedeberg]kohl2016directed
authorM. Kohl,
  authorR. Capellmann,
  authorM. Laurati,
  authorS. Egelhaaf, and
  authorM. Schmiedeberg,
  journalNat. Commun. volume7
  (year2016).

[Verhaegh et al.(1997)Verhaegh, Asnaghi,
  Lekkerkerker, Giglio, and Cipelletti]verhaegh1997
authorN. A. M. Verhaegh,
  authorD. Asnaghi,
  authorH. N. W. Lekkerkerker,
  authorM. Giglio, and
  authorL. Cipelletti,
  journalPhysica A volume242,
  pages104 (year1997).

[Pham et al.(2002)Pham, Puertas,
  Bergenholtz, Egelhaaf, Moussaıd, Pusey, Schofield, Cates, Fuchs, and
  Poon]pham2002multiple
authorK. N. Pham,
  authorA. M. Puertas,
  authorJ. Bergenholtz,
  authorS. U. Egelhaaf,
  authorA. Moussaıd,
  authorP. N. Pusey,
  authorA. B. Schofield,
  authorM. E. Cates,
  authorM. Fuchs, and
  authorW. C. Poon,
  journalScience volume296,
  pages104 (year2002).

[Baxter(1968)]baxter1968percus
authorR. Baxter, journalJ.
  Chem. Phys volume49, pages2770
  (year1968).

[Buzzaccaro et al.(2007)Buzzaccaro,
  Rusconi, and Piazza]buzzaccaro2007sticky
authorS. Buzzaccaro,
  authorR. Rusconi, and
  authorR. Piazza,
  journalPhys. Rev. Lett. volume99,
  pages098301 (year2007).

[He and Thorpe(1985)]he1985elastic
authorH. He and
  authorM. Thorpe,
  journalPhys. Rev. Lett. volume54,
  pages2107 (year1985).

[Hinrichsen(2000)]hinrichsen2000non
authorH. Hinrichsen,
  journalAdv. Phys volume49,
  pages815 (year2000).

[Bodnár et al.(1996)Bodnár,
  Dhont, and Lekerkerker]bodnar1996
authorI. Bodnár,
  authorJ. K. G. Dhont,
  and authorH. N. W.
  Lekerkerker, journalJ. Phys. Chem.
  volume100, pages19614
  (year1996).

[Royall et al.(2007)Royall, Louis, and
  Tanaka]royall2007measuring
authorC. P. Royall,
  authorA. A. Louis, and
  authorH. Tanaka,
  journalJ. Chem. Phys volume127,
  pages044507 (year2007).

[Lekkerkerker et al.(1992)Lekkerkerker,
  Poon, Pusey, Stroobants, and Warren]lekkerkerker1992
authorH. N. W. Lekkerkerker,
  authorW. C. K. Poon,
  authorP. N. Pusey,
  authorA. Stroobants,
  and authorP. B.
  Warren, journalEurophys. Lett.
  volume20, pages559 (year1992).

[Zhang et al.(2013)Zhang, Royall, Faers,
  and Bartlett]zhang2013
authorI. Zhang,
  authorC. P. Royall,
  authorM. A. Faers, and
  authorP. Bartlett,
  journalSoft Matter volume9,
  pages2076 (year2013).

[Leocmach and Tanaka(2013)]leocmach2013novel
authorM. Leocmach and
  authorH. Tanaka,
  journalSoft Matter volume9,
  pages1447 (year2013).

[Poon et al.(2012)Poon, Weeks, and
  Royall]poon2012measuring
authorW. C. Poon,
  authorE. R. Weeks, and
  authorC. P. Royall,
  journalSoft Matter volume8,
  pages21 (year2012).

[Royall et al.(2013)Royall, Poon, and
  Weeks]royall2013search
authorC. P. Royall,
  authorW. C. Poon, and
  authorE. R. Weeks,
  journalSoft Matter volume9,
  pages17 (year2013).

[Turci and Royall(2016)]turci2016crystallisation
authorF. Turci and
  authorC. P. Royall,
  journalJ. Stat. Mech. Theory Exp volume2016,
  pages084004 (year2016).

[Brostow et al.(1978)Brostow, Dussault,
  and Fox]brostow1978construction
authorW. Brostow,
  authorJ.-P. Dussault,
  and authorB. L. Fox,
  journalJ. Comput. Phys volume29,
  pages81 (year1978).

[Malins et al.(2013)Malins, Williams,
  Eggers, and Royall]malins2013identification
authorA. Malins,
  authorS. R. Williams,
  authorJ. Eggers, and
  authorC. P. Royall,
  journalJ. Chem. Phys volume139,
  pages234506 (year2013).

[Binder(1981)]binder1981finite
authorK. Binder,
  journalZeitschrift für Physik B Condensed Matter
  volume43, pages119 (year1981).

[Watanabe et al.(2012)Watanabe, Ito, and
  Hu]watanabe2012phase
authorH. Watanabe,
  authorN. Ito, and
  authorC.-K. Hu, journalJ.
  Chem. Phys volume136, pages204102
  (year2012).

[Rovere et al.(1993)Rovere, Nielaba, and
  Binder]rovere1993simulation
authorM. Rovere,
  authorP. Nielaba, and
  authorK. Binder,
  journalZeitschrift für Physik B Condensed Matter
  volume90, pages215 (year1993).

[Sengupta et al.(2000)Sengupta, Nielaba,
  Rao, and Binder]sengupta2000elastic
authorS. Sengupta,
  authorP. Nielaba,
  authorM. Rao, and
  authorK. Binder,
  journalPhys. Rev. E volume61,
  pages1072 (year2000).

[Galkin and Vekilov(2000)]galkin2000control
authorO. Galkin and
  authorP. G. Vekilov,
  journalProceedings of the National Academy of Sciences
  volume97, pages6277 (year2000).

[Fortini et al.(2008)Fortini, Sanz, and
  Dijkstra]fortini2008crystallization
authorA. Fortini,
  authorE. Sanz, and
  authorM. Dijkstra,
  journalPhys. Rev. E volume78,
  pages041402 (year2008).

[Savage and Dinsmore(2009)]savage2009
authorJ. R. Savage and
  authorA. D. Dinsmore,
  journalPhys. Rev. Lett. volume102,
  pages198302 (year2009).

[Ohtsuka et al.(2008)Ohtsuka, Royall, and
  Tanaka]ohtsuka2008
authorT. Ohtsuka,
  authorC. P. Royall,
  and authorH. Tanaka,
  journalEurophys. Lett. volume84,
  pages46002 (year2008).

[Rice et al.(2012)Rice, Roth, and
  Royall]rice2012
authorR. Rice,
  authorR. Roth, and
  authorC. P. Royall,
  journalSoft Matter volume8,
  pages1163 (year2012).

[Varga and Swan(2016)]varga2016hydrodynamic
authorZ. Varga and
  authorJ. Swan,
  journalSoft matter volume12,
  pages7670 (year2016).

[Doye et al.(1995)Doye, Wales, and
  Berry]doye1995effect
authorJ. P. Doye,
  authorD. J. Wales, and
  authorR. S. Berry,
  journalJ. Chem. Phys volume103,
  pages4234 (year1995).

[Malins et al.(2009)Malins, Williams,
  Eggers, Tanaka, and Royall]malins2009
authorA. Malins,
  authorS. R. Williams,
  authorJ. Eggers,
  authorH. Tanaka, and
  authorC. P. Royall,
  journalJ. Phys.: Condens. Matter volume21,
  pages425103 (year2009).

[Karayiannis et al.(2011)Karayiannis,
  Malshe, de Pablo, and Laso]karayiannis2011fivefold
authorN. C. Karayiannis,
  authorR. Malshe,
  authorJ. J. de Pablo,
  and authorM. Laso,
  journalPhys. Rev. E volume83,
  pages061505 (year2011).

[Taffs et al.(2013)Taffs, Williams,
  Tanaka, and Royall]taffs2013structure
authorJ. Taffs,
  authorS. R. Williams,
  authorH. Tanaka, and
  authorC. P. Royall,
  journalSoft Matter volume9,
  pages297 (year2013).

[Royall et al.(2015b)Royall,
  Eggers, Furukawa, and Tanaka]royall2015
authorC. P. Royall,
  authorJ. Eggers,
  authorA. Furukawa, and
  authorH. Tanaka,
  journalPhys. Rev. Lett. volume114,
  pages258302 (year2015b).



