DOI:
10.1039/D3SM01314F
(Paper)
Soft Matter, 2023,
19, 9531-9540
Glassy phases of the Gaussian core model
Received
2nd October 2023
, Accepted 26th November 2023
First published on 28th November 2023
Abstract
We present results from molecular dynamics simulations exploring the supercooled dynamics of the Gaussian Core Model in the low- and intermediate-density regimes. In particular, we analyse the transition from the low-density hard-sphere-like glassy dynamics to the high-density one. The dynamics at low densities is well described by the caging mechanism, giving rise to intermittent dynamics. At high densities, the particles undergo a more continuous motion in which the concept of cage loses its meaning. We elaborate on the idea that these different supercooled dynamics are in fact the precursors of two different glass states.
I. Introduction
Soft colloids encompass a large variety of nano- or mesoscopic aggregates undergoing thermal motion in a solvent. In most cases, these are polymer-based assemblies of various architectures and connectivities, such as linear chains, stars and rings, dendrimers, block-copolymer micelles, but also cross-linked nano- or microgels, which feature a high degree of flexibility and deformability.1 The degree of softness is conveniently expressed by the ratio Eel/kBT of their elastic deformation energy upon small overlaps or indentations, Eel, over the thermal energy kBT, which can span several orders of magnitude.1,2 Their suspensions exhibit structural and dynamical anomalies,3–6 which accompany rich types of thermodynamic behavior, such as reentrant melting,7–11 and clustering.12,13 The deformability of these soft colloids manifests itself under flow as particle elongation, tumbling and shear-thinning.14–18 Akin to hard colloids, in the supercooled regime, soft colloids exhibit a sharp increase of the equilibrium relaxation time and heterogeneous dynamics.1,19–21 However, at very high concentrations, soft colloids can reach a particular aging regime characterized by an intermittent release of internal stresses which coincides with the onset of an anomalous decrease in local order.22 Colloidal softness has been also related to the microscopic origin leading to the validity of the Stokes–Einstein (SE) relation for degrees of metastability for which it normally breaks down in the case of hard colloidal and molecular systems.23
From a theoretical point of view, the effective interaction v(r) between a pair of such colloids does not grow fast as r → 0, so that the integral
is finite. In some cases, these are effective interactions between the centers of mass of open, fractal, fully penetrable objects, for which v(r) remains finite (free of divergence) even if the separation r vanishes. Denoting with
(k) the Fourier transform of v(r), such effective interactions are referred to as Q+ potentials if
(k) is positive definite and as Q± potential if
(k) attains both positive and negative values, typically in an oscillatory fashion.24 Systems belonging to the Q+ class undergo a reentrant fluid–crystal–fluid transition at low temperature and high density and they possess a maximum freezing temperature, beyond which no crystallization is possible whereas systems of the Q± class bond together forming cluster crystals, where each lattice site is occupied by several overlapping particles.13,24
Within the realm of Q+ potentials a prominent role is played by the Gaussian Core Model (GCM) introduced by Stillinger in the 70's.7 The GCM is one of the simplest models for the description of systems such as polymer or dendrimer solutions25,26 and it entails an inter-particle Gaussian-shaped potential
| v(r) = ε exp[−(r/σ)2], | (1) |
where
ε and
σ are the energy and length scales. In what follows, we express time in units of

, where
m is the mass of each particle. We further define a reduced density
ρ →
ρσ3 and a reduced temperature
T →
kBT/
ε, where
ρ =
N/
V is the number density,
kB is Boltzmann's constant, and
T the absolute temperature. Whereas at low temperatures and low densities, the equilibrium properties of the GCM can be described by an effective hard-sphere mapping,
7 at high densities a mean-field description sets in, giving rise to re-entrant melting below a threshold upper freezing temperature
Tu = 8.74 × 10
−3, above which it remains fluid at all densities.
7,27,28 The equilibrium phase diagram of the GCM, along with two predictions for its vitrification line, arising from two different approximations, is shown in
Fig. 1. Another prominent member of the
Q+ -class is the Hertzian potential,
v(
r) =
ε(1 −
r/
a)
5/2Θ(1 −
r/
a), which models the effective interaction between elastic spheres of diameter
a.
8 The general features of the phase diagram of the Hertzian spheres are similar to those of the GCM. It must be emphasized, however, that the former has finite support, as expressed by the Heaviside function
Θ(1 −
r/
a), whereas the GCM is non-vanishing for arbitrarily large values of
r. In this work, we focus our attention to the glassy states of the GCM and we establish that as those are approached from the supercooled liquid, there are two distinct types of single-particle dynamics in the latter: at low densities, intermittent dynamics akin to the standard caging scenario is present whereas at intermediate densities the dynamics becomes more continuous and the classical cage picture is gradually lost as density grows.
 |
| Fig. 1 Phase diagram of the GCM. The black solidification line and the red glass line from the Replica Theory (RT) are taken from the literature, in particular from ref. 28 and 33, respectively. The green diamonds indicate the state points with the same diffusion coefficient, D ≃ 1.4 × 10−4; the green dashed line is a guide to the eye. The blue glass line is calculated from Mode Coupling Theory (MCT) as described in Materials and Methods. The MCT line is recovered at densities ρ ≳ 1.00, where it follows a monotonically decreasing trend with density;29,30 see main text. | |
II. Problem setup
At low densities, the GCM can be mapped onto an effective hard sphere potential,7 and its crystallization properties can be understood in such a context. The high-density part of the phase behavior of the GCM is more challenging. Previous studies29,30 have shown that at high densities nucleation is strongly suppressed and that Mode-Coupling Theory (MCT) provides an accurate description of the structural arrest of the GCM into an amorphous, glassy state. Moreover, it was found that the glass state at high densities displays strong dynamic fluctuations and a nearly Gaussian distribution of single particle displacements, features compatible with a geometric transition.31 Such transition refers to a change in the topology of the rugged free energy landscape; in particular, below the transition temperature, the potential energy barriers become much larger than the available thermal energy, and the system is trapped close to local minima of the energy landscape.32
The low- and intermediate-density glassy regimes of the GCM remain largely unexplored. Assuming that the low-density vitrification scenario follows the hard-sphere paradigm, it is then particularly interesting to ask the question as to how the vitrification scenario evolves towards the high-density regime and whether dynamically distinct glassy states exist in different density regimes of the GCM. A recent theoretical study33 showed an unexpected density dependence of the glassy behavior of GCM particles, see Replica Theory (RT) line in Fig. 1. Similarly to the equilibrium crystallization behavior, the glass line shows a re-entrance upon increasing the density. However, at moderate densities, the characteristic order parameter at constant density displays sudden jumps when increasing the temperature. This trend suggests a transition between two different glasses, a continuous and a discretized one. In particular, the emergence of a discretized glass has been associated with the formation of out-of-equilibrium local aggregates. Indeed, as mentioned above, the GCM is a Q+ potential for which no clusters form at equilibrium. In contrast, for ultrasoft particles belonging to the Q± class for which cluster formation is an equilibrium phenomenon, the emergence of cluster-glasses has been recently reported.34–37 We note that in Fig. 1, the MCT-vitrification line stops at ρ ≳ 0.40 and at about the same density the RT-vitrification line shows non-monotonic behavior with density. The reason for the former is a convergence loss of the Ornstein–Zernike equation in the hypernetted chain (HNC) approximation for the one-component GCM, which is however recovered at much higher densities, ρ ≳ 1.00. The RT, being a two-component approach, does converge in this region and results in the aforementioned non-monotonic behavior. Whereas it is an open question whether this behavior is connected with the convergence problems of the HNC in that region, the discretized glass predicted by the RT occurs already at lower densities, and thus the question of whether a distinct arrested state exists there is independent of the HNC-convergence issues.
The goal of this work is to characterize the transition from low- to high-density glass from a dynamic point of view, focusing on the study of the supercooled regime. Indeed, when approaching the glass transition the system enters into a supercooled regime which represents in all respects a precursor of the glass.38 Thus, we expect to observe differences between the two states already at the level of supercooled (glassy) dynamics. The supercooled regime of canonical glass-formers is usually described in terms of the caging mechanism: each particle experiences trapping due to the neighboring particles that effectively create a cage around it; eventually, the fluctuations allow the particle to escape this local cage and move to the next one. The lower the temperature the more difficult it will be for the particle to escape from the cage. In these terms, the glass transition can be thought of as a localization transition. The cage size is related to the average inter-particle distance, which in turn depends on the density of the system. Such a mechanism is accurate for systems characterized by a harshly repulsive inter-particle potential, such as hard-sphere or Lennard-Jones systems.39,40 However, when dealing with bounded potentials, that is with potentials that do not diverge when two particles are at full overlap, and in the presence of long tails, this mechanism can break down. In particular, we show that for the GCM at intermediate densities, the idea of slow dynamics based on the concept of caging must be revised.
III. Results
Following the phase diagram in Fig. 1, we simulate the glassy dynamics of the GCM at different densities. As mentioned above, at high densities the one-component GCM vitrifies and thus, it is possible to approach the supercooled regime directly, without running into crystallization issues. This is not the case for low and intermediate densities, for which the one-component system would crystallize upon cooling. Therefore, in our simulations, we follow a random pinning procedure41–44 and freeze a fraction f of the particles to avoid crystallization and be able to approach the deeply supercooled regime. The latter is accomplished by considering that increasing the pinning fraction has the same effects as decreasing the temperature for the bulk system.45 All results reported below are calculated for f = 10%, taking into account the mobile particles only, and were averaged over at least three different realizations of the pinning disorder (see Materials and Methods for more details on the simulation protocol).
We focus our analysis on four different densities, i.e., ρ = 0.10, 0.15, 0.40, and 1.00, in order to characterize the transition between low- and high-density glassy dynamics. For each density, we spanned over various different temperatures and selected those at which all four systems display the same diffusivity at long times, as can be seen at the long-time behavior of the particle displacements shown in Fig. 2(a). We select points that lie along such an isodiffusivity line, for three different reasons. First, it has already been established for the GCM that the shape of these lines is very similar to that of lines on which the main peak of the structure factor or of the radial distribution function is constant, signifying very similar pair correlations for all state points lying on them.46,47 Second, not only the peak heights themselves but the whole structure factors of points on the isodiffusivity lines can be essentially rescaled on a single curve, as shown in Fig. 2(b) and further discussed below. Finally, there is a strong correlation between the isodiffusivity and the vitrification lines,19 the former acting as a precursor to the latter and thus any two points on an isodiffusivity lines being ‘equidistant’ from the glass line in terms of the temperature quench needed to vitrify the system.
 |
| Fig. 2 (a) Main plot: the mean squared displacement 〈Δr2(t)〉 of the GCM as a function of time, calculated at four different isodiffusive points, as indicated in the legend. Inset: Corresponding local exponent of the MSD as defined in (6). (b) structure factors S(q) calculated at the same isodiffusive points as panel (a); in the horizontal axis, the wavenumbers q of each structure factor are rescaled over the value related to the corresponding first peak (qpeak). | |
To characterise particle motion, we first calculated the Mean-Square Displacements (MSD) as
|  | (2) |
where
Nm is the number of mobile particles and the brackets 〈⋯〉 denote an average over all particles, with initial positions
ri(
t0) at time
t0 and
ri(
t +
t0) at a time interval
t later. We also average over different initial times
t0.
The selected temperature values are reported in Fig. 2(a), in which the MSD is plotted, clearly showing the same long-time diffusivity for the four systems. We have also calculated the equal-time static structure factor
|  | (3) |
which is shown in
Fig. 2(b). After suitable rescaling due to the different densities, the overall structural properties of the four systems are comparable. Pairs of states that lie on the same isodiffusivity line for a given
T-value can thus have their pair correlation functions be mapped onto one another upon rescaling the lengths by the average interparticle distance
a; we term such states conjuguate pairs. This is a kind of extension into finite temperatures and into the fluid state of the exact duality relation at
T = 0 that connects the ground-state energies of the fcc- and bcc-lattices,
Ifcc(
a) and
Ibcc(
a), which reads as:
48 |  | (4) |
This property, together with the same long-time diffusivity, would suggest that a re-mapping of all systems is possible and that the particle dynamics of the associated systems can also be collapsed onto a single curve upon suitable rescalings. This, however, is not the case, as we discuss below.
We consider next the intermediate scattering function (ISF), defined as
|  | (5) |
which is reported in
Fig. 3, revealing important differences in the relaxation dynamics. In particular, at the two lower densities, a clear two-step relaxation can be discerned, which becomes smoothed out at the intermediate one and practically disappears at the highest of the four, indicating the lack of a caging mechanism for the latter. Therefore, the ISF displays a faster relaxation for high density, while for low and intermediate densities the two-step behavior emerges, leading to an effectively slower relaxation. Due to the presence of the pinned particles, it is expected that the coherent part of the ISF does not decay completely to zero at long times.
44,49 Before this expected long-time behavior is reached, both the self and the coherent parts of the ISF display a similar trend, with no indication of decoupling between them, at least in the regime investigated within this work. There is, therefore, one single relaxation time associated with both the incoherent and the coherent intermediate scattering functions and not two separate ones, as is the case with other single-component ultrasoft systems, such as semiflexible minirings.
50 We obtained the relaxation time
τα by setting
Fs(
q,
t =
τα) = 1/
e and compared it with the diffusion coefficient extracted from the MSD in
Fig. 2(a) in order to study the violation of the SE relation (
Dτα ∼ const). It is known that the SE relation breakdown is very weak at high density;
30 the inset reported in
Fig. 3 shows that this is not the case for low density (
ρ = 0.15), at which we observe strong deviations from the SE relation. At intermediate density (
ρ = 0.40) the SE breakdown becomes more moderate, approaching the well known results of weak violations at high density.
30
 |
| Fig. 3 The self-intermediate scattering function Fs(qpeak,t) for the four different isodiffusive points, evaluated at the corresponding wavenumber qpeak where the structure factor of each system has its main peak. The inset highlights the SE relation breakdown at low and intermediate density. The quantity Dτα is normalized over its value at high temperature, indicated as (Dτα)ref. | |
The differences between low- and intermediate densities emerging from the the study of the ISF can in fact already be noticed in Fig. 2(a), focusing on the behavior of the MSD at times between the ballistic and diffusive regimes. For the lower densities, a much stronger caging effect, indicated by the development of a plateau in the MSD, emerges than for the higher ones. It is worth mentioning that the lack of a plateau at the higher density precisely compensates the fact that the ballistic motion is faster for the lower densities, so that when the former can enter a plateau while the motion at the higher density catches up. Consequently, all MSD's enter their diffusive regime at the same distance squared and at the same time, following thereafter the same diffusive pattern. These differences become quantitative by considering the local exponent
of the MSD, defined as
|  | (6) |
and shown in the inset of
Fig. 2(a). Whereas the local exponents for the two lowest densities show a clear transition from ballistic

to a diffusive

regimes through an intermediate plateau

, both the ballistic regime and the clear intermediate plateau disappear at the high-density part of the isodiffusivity line. This is a clear indication that at short-to-intermediate scales, the two motions differ, a prediction to be quantified and analyzed below by statistically analyzing individual particle trajectories.
In Fig. 4 we report typical single-particle displacements,
. By looking at the left panels of Fig. 4 it is clear that the dynamics shifts from an intermittent-like behavior at low densities to a more continuous one at high densities. Different statistical methods have been used in literature to detect jumps between cages and thus characterize the intermittent dynamics.40,51–54 In general, most of these methods consider a designated estimator for the distance covered between two points in time and compare it with a specific cut-off distance related to the cage size. Such methods, which are defined in a differential form, are rather sensitive to noise. Moreover, they rely on a clear estimate of the cage size which, in the case of soft systems, is not always easy to define. To overcome these issues, we make use of a recently developed analysis based on the Local Convex Hull (LCH) method, which turns out to be much less sensitive to fluctuations and does not require an a priori knowledge on the cut-off distance.55
 |
| Fig. 4 Typical single-particle displacements (left column) and corresponding SV(t) time series from the LCH analysis (right column) calculated as described in Materials and Methods. The black solid line indicates the threshold dividing slow and fast phases and is calculated as the average over the the whole time series SV(t), that is . | |
The main idea behind this analysis is to use geometric properties of the smallest convex shape (precisely the LCH) enclosing a small set of trajectory points to estimate the space explored by each particle in a specific time window, as shown in Fig. 5. More specifically, our analysis focuses on the study of the LCH volume SV(t) which, with respect to other geometric quantities such as the diameter, is more sensitive to changes in the dimensionality and anisotropy of the particle motion.55 In Fig. 4, together with the single-particle displacements, we report also the corresponding time series SV(t) calculated from the LCH method as described in Materials and Methods. We can observe that, if the particle motion has an intermittent-like behaviour, SV(t) will display few and high peaks, while, if the particle motion follows a more continuous trend, SV(t) will mostly oscillate around its single-particle average value
with multiple, lower peaks. Such trend suggests that performing a statistical analysis of the SV(t) peaks can help in classifying different dynamical behaviors. Then, for each time series we find the peak locations and corresponding peak values
. Moreover, we identify with ΔtSP the duration of the so called slow phases, that is the time during which SV(t) stays continuously below the threshold
before crossing it. Then, each trajectory will display a certain number of slow phases and peak values.
 |
| Fig. 5 Single-particle trajectory and corresponding LCHs calculated for two sets of points centered in different time instants (red points). The volume of the LCH in (a) is clearly smaller than the one in (b), suggesting that a jump between two local cages can be identified as a peak in the time series of SV(t) as it is visible in Fig. 4. | |
In Fig. 6 we report the probability distributions of (i) ΔtSP, (ii) the peak height evaluated from the threshold value (e.g.
), and (iii) the number of peaks. We can see immediately that p(ΔtSP) classifies the systems into two different dynamics: the distribution presents a fatter tail for the two lower density systems with respect to the higher-density systems. In addition, by looking at p(
) we observe that the peak height is more likely to assume larger values for the lower density systems than for the higher density ones. However, in this case, there is not a full rescale of the high-density systems, suggesting that ρ = 0.40 still belongs to a transition phase between the two dynamical regimes. We emphasize that the quantity p(
) has been rescaled by the average value
in order to eliminate trivial contributions stemming from the density-dependence of the volume explored by each particle (see also Fig. 7). Finally, the distribution of the number of peaks
complements the information provided by p(ΔtSP) suggesting an intermittent-like motion when fewer peaks (and longer slow phases) are detected and a more continuous one when a larger number of peaks (and shorter slow phases) is observed. Indeed,
shows that the system with the highest density is shifted towards larger values of peak numbers with respect to the two systems with lower densities, which follow a similar distribution centered around smaller values; once again, the system with ρ = 0.40 displays an intermediate behavior confirming that at this density the system is in a transition phase between the two dynamical regimes.
 |
| Fig. 6 Probability distribution of (a) slow phase duration ΔtSP, (b) peak height with respect to the threshold value and (c) number of peaks. The distribution in (b) is rescaled by the average value which clearly depends on the density (see Fig. 7). | |
 |
| Fig. 7 Probability distribution of obtained by considering all the values extracted from each single-particle trajectory as indicated in Materials and methods. | |
We further investigate the particle-to-particle variability of the threshold value
, which indicates the average volume explored by each particle within the simulation time window. To do so, we extract the value
for each (mobile) particle and then build the histogram, as reported in Fig. 7. It can be seen that the distribution
looks quite narrow for the highest density, suggesting a more homogeneous dynamics. Conversely,
for the lower density is much broader, implying the presence of slow and fast particles in agreement with the concept of dynamic heterogeneity typical of canonical supercooled liquids. Of particular importance is the fact that for the lower densities the distribution has contributions even at
, suggesting that, within our time window, there are particles that do not move at all (or only very little) from their initial cage. This is not the case for the higher densities, for which the distribution vanishes for small values of
.
Additional corroboration for the gradual disappearing of the standard cage-escape mechanism for the supercooled GCM-liquid at intermediate densities is offered by considering the single-particle displacement distribution P(Δr,t) and the corresponding non-Gaussian parameter α2(t) = 3〈Δr4(t)〉/5[〈Δr2(t)〉]2 − 1. The latter is used to quantify the deviation of particles displacements from a Gaussian distribution, identified by α2(t) = 0. When dealing with supercooled liquids, such deviations of this quantity from zero are usually attributed to the presence of dynamic heterogeneity in the system. However, in the context of the high-density GCM, it has been shown that small values of α2(t) are compatible with strong dynamical heterogeneities emerging from a mean-field, geometric glass transition.31 In Fig. 9 we show the calculated non-Gaussian parameter for the four isodiffusivity state points, finding that those of the low-density points differ drastically from those of the high-density points. In particular, there is a non-monotonic behavior of the curves and of their maximum values, which occurs roughly at the end of the caging time for ρ = 0.10 and ρ = 0.15 and somewhat earlier for ρ = 0.40 and ρ = 1.00. Thus, α2(t) follows the same trend in density as the isodiffusivity line and other quantities characterizing the system, see Fig. 1 and ref. 4, 46 and 56 Single-particle motions tend thus to become more and more Gaussian as the density grows, in agreement with the absence of a bimodal (cage/cage escape jump) distribution of the self-van Hove function31 and thus with the gradual disappearance of the cage-hopping dynamics characteristic of the low-density supercooled GCM fluid. In fact, we explicitly confirm the suppression of hopping dynamics at intermediate and higher densities along our isodiffusivity line in Fig. 8, demonstrating indeed that the standard caging mechanism of dynamic slowing down in the supercooled liquid is valid only on the low-density side of the vitrification line.
 |
| Fig. 8 Evolution in time of the single-particle displacement distributions for the four isodiffusivity state points: (a) t = 102, (b) t = 103, (c) t = 104, (d) t = 105. In particular, we calculated the probability distributions of the logarithm of single-particle displacements which allow us to highlight the hopping motion, when present. Indeed, in panel (c) we can clearly see the emergence of a bimodal distribution for all but the highest density system. | |
 |
| Fig. 9 The non-Gaussian parameter α2(t) of the GCM-supercooled fluids at the four isodiffusivity state points considered. | |
Finally, in Fig. 10 we report the non-ergodicity factor calculated within the MCT framework. As T decreases, the low-density curves tend to behave similarly and the same happens for the high-density graphs. In addition, the non-ergodocity factors for q/qpeak > 1 are also very similar to one another for all the points located along the isodiffusivity lines. Significant differences appear, on the other hand, for the nonergodocity factors at large values of the wavelength. Here, the intermediate-density curves feature much lower values of the ϕ(q) than their low-density counterparts. Our results are in this way precursors of the findings by previous works at the high-density part of the GCM, referring to the mean-field glass for which the MCT becomes vey accurate.29–31 There, and for density ρ = 2, it was found that the low-vavenumber values of the nonergodicity factor become extremely small and that a decoupling between collective- and self-dynamics at these scales follows, as a result of the extremely low compressibility of the GCM under such conditions. The spectacular suppression of the nonergodicity plateau at small q-values for the intermediate-density GCM glass expresses an anticipation of the ideal MCT-glass at high densities, which features long-wavelength mobility modulations.29–31
 |
| Fig. 10 Non-ergodicity factor ϕ(q) obtained from MCT (see materials and methods). We compare same temperatures at low density (orange curves) and intermediate density (red curves). | |
IV. Discussion and conclusions
The intermittent dynamics identified for the low-density GCM in the supercooled regime confirms the picture provided by the hard-sphere re-mapping, which appears to remain valid also for the glassy regime. Indeed, in this regime, the GCM displays a glassy behavior that is in agreement with the one of canonical HS-like glass formers. In particular, the dynamics is fully described by the hopping between local cages leading to standard features such as the emergence of the plateau in the MSD, dynamical heterogeneity, and non-Gaussianity. All these features have been re-classified in this work by making use of the LCH analysis. Upon increasing the density, the GCM undergoes a transition towards a different glassy state. The intermittent-like behavior is taken over by a more continuous more Gaussian dynamics. Our analysis shows how the slowing down at high densities of the GCM cannot be explained using the standard caging mechanism. As a matter of fact, the higher the density the less confining the potential becomes due to the increasing contribution of neighbouring particles to the energy landscape. The dynamics in this regime approaches more and more a mean-field-like picture, as already anticipated in.31 With our single-particle analysis, we were able to fully capture the transition between these two glassy states (HS-like at low densities and mean-field-like at high densities) and to show how the intermediate density regime is characterized by a smooth transition between the two dynamics. Whether this transition in the supercooled fluid regime is the echo of an underlying glass-glass transition between two arrested states, akin to the distinct glassy states encountered, e.g., in colloid-polymer mixtures,57 is an open problem and it will be the subject of further investigations.
V. Materials and methods
Molecular dynamics simulations
Molecular dynamics (MD) simulations were performed using the open-source package LAMMPS58 for a cubic box with periodic boundary conditions and several combinations of density ρ = N/V and temperature T, where N and V are the total number of Gaussian particles and the volume of the box, respectively. The equations of motion were integrated using the velocity Verlet scheme,59 with a time step Δt/τ = 0.01, where
and m the mass unit for the particles.
Due to the propensity of the system to crystallize at low temperatures, it is required to add some frustration degree into the system, which helps particles to remain in a disordered phase as temperature decreases. In this work, frustration is introduced through random pinning that avoids the inclusion of randomness in the interaction potential or compositional disorder.41–45
The initial configuration of the system was generated in a multiple steps process: at first, Np = fN particles were randomly placed in the simulation box and an equilibration process took place at a relatively high temperature (T = 0.01). Then they were permanently pinned and the remaining Nm = N − Np particles (the mobile ones) were randomly inserted in the box preventing excessive overlapping with the previously inserted particles. By setting the target temperature T, an equilibration run was performed for the mobile particles in which the system was thermalized by means of the Nose–Hoover thermostat. To consider the effect of the pinning protocol, a second pinning method was also employed. In that case, the whole system is equilibrated at a high temperature, then the fraction f of particles is permanently frozen and the system is lastly quenched to the target temperature for the equilibration of the mobile particles to occur.
We worked with a total number of particles N = 3456 and f = 0.10. The typical equilibration time for the mobile particles ranged from 2 × 106 to 107 time steps. After reaching a steady state, indicated by the absence of any drift in internal energy and pressure, a production run was performed, which typically ranged from 107 to 108 time steps. Provided that all average measurements for given density and temperature are performed over both thermal fluctuations and different realizations of the pinning disorder, it is expected that the thermodynamics of the mobile particles remains unperturbed.42
Local convex hull analysis
To differentiate between different dynamics phases we analyse the single-particle trajectories making use of the Local Convex Hull (LCH) method developed in ref. 55. The convex hull of a finite set of points
is the minimal convex shape that encloses all the points x1,…,xn. We consider the volume QV(i) of the LCH computed over 2τ0 + 1 trajectory points centered around xi. Then, we classify the point xi by using all the estimators to which it contributes (in total 2τ0 + 1), that is by introducing the discriminator: |  | (7) |
The LCH method consists in: (1) mapping each trajectory into a one-dimensional time series by computing the discriminator SV(n) for all points xn; and (2) selecting a threshold SV such that the points xi with SV(i) > SV are classified into the fast phase whereas the points xi with SV(i) ≤ SV into the slow phase. Following ref. 55 we set the threshold SV to be the average value of SV(n) over the single trajectory, namely |  | (8) |
The parameter τ0 was chosen empirically by looking at the MSD curves and selecting the time at which the system escapes from the plateau, generally used as an estimate of the caging time. For our system in particular we have τ0 = 103 (in simulation time units). To calculate the LCH, we used the available algorithm in Python based on the Qhull library.
Mode coupling theory
To quantify the liquid-to-glass transition of the GCM system, the non-ergodicity factor ϕ(q) was evaluated, which is defined as the long-time limit of the density autocorrelation function, i.e., |  | (9) |
where
and the sum is performed over the coordinates rj(t) of all particles in the system.
In the case ϕ(q) ≠ 0, the system is considered non-ergodic and its state is identified as glassy, whereas ϕ(q) = 0 corresponds to an ergodic fluid. Given the structural data obtained by solving the Ornstein–Zernike equation through the hypernetted-chain closure (HNC), the calculation of the non-ergodicity factor is readily achieved within the framework of the Mode Coupling Theory (MCT). According to it, ϕ(q) fulfills the self-consistent equation:60
|  | (10) |
where the kernel

can be expressed entirely in terms of the Fourier transform of the direct correlation function
ĉ(
k),
i.e.,
|  | (11) |
and
S(
k) = [1 −
ρĉ(
k)]
−1.
Conflicts of interest
There are no conflicts to declare.
Acknowledgements
We thank S. Prestipino for providing the data from ref. 28 shown in Fig. 1 and E. Zaccarelli for helpful discussions. V. S. acknowledges support of the European Commission through the Marie Skłodowska-Curie COFUND project REWIRE, Grant Agreement No. 847693. M. C. thanks Universidad Antonio Nariño for financial support through Project No. 2022202.
References
- D. Vlassopoulos and M. Cloître, Tunable rheology of dense soft deformable colloids, Curr. Opin. Colloid Interface Sci., 2014, 19, 561–574 CrossRef CAS.
- J. Riest, L. Athanassopoulou, S. A. Egorov, C. N. Likos and P. Ziherl, Elasticity of polymeric nanocolloidal particles, Sci. Rep., 2015, 5, 15854 CrossRef PubMed.
- M. Watzlawek, H. Löwen and C. N. Likos, The anomalous structure factor of dense star polymer solutions, J. Phys.: Condens. Matter, 1998, 10, 8189–8205 CrossRef CAS.
- W. P. Krekelberg, T. Kumar, J. Mittal, J. R. Errington and T. M. Truskett, Amonalous structure and dynamics of the Gaussian-core fluid, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2009, 79, 031203 CrossRef PubMed.
- M. J. Pond, J. R. Errington and T. M. Truskett, Communication: Generalizing Rosenfeld's excess-entropy scaling to predict long- time diffusivity in dense fluids of Brownian particles: From hard to ultrasoft interactions, J. Chem. Phys., 2011, 134, 081101 CrossRef PubMed.
- M. J. Pond, J. R. Errington and T. M. Truskett, Mapping between long-time molecular and Brownian dynamics, Soft Matter, 2011, 7, 9859–9862 RSC.
- F. H. Stillinger, Phase transitions in the Gaussian core system, J. Chem. Phys., 1976, 65, 3968–3974 CrossRef CAS.
- J. C. Pamies, A. Cacciuto and D. Frenkel, Phase diagram of Hertzian spheres, J. Chem. Phys., 2009, 131, 044514 CrossRef PubMed.
- L. Berthier, A. J. Moreno and G. Szamel, Increasing the density melts ultrasoft colloidal glasses, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2010, 82, 060501 CrossRef PubMed.
- N. Gnan and E. Zaccarelli, The microscopic role of deformation in the dynamics of soft colloids, Nat. Phys., 2019, 15, 683 Search PubMed.
- D. Parisi, D. Truzzolillo, A. H. Slim, P. Dieudonné-George, S. Narayanan and J. C. Conrad,
et al., Gelation and Re-entrance in Mixtures of Soft Colloids and Linear Polymers of Equal Size, Macromolecules, 2023, 56, 1818–1827 CrossRef CAS PubMed.
- D. Lenz, R. Blaak, C. N. Likos and B. M. Mladek, Microscopically resolved simulations prove the existence of soft cluster crystals, Phys. Rev. Lett., 2012, 109, 228301 CrossRef PubMed.
- E. Stiakakis, N. Jung, N. Adžić, T. Balandin, E. Kentzinger and U. Rücker,
et al., Self assembling cluster crystals from DNA based dendritic nanostructures, Nat. Commun., 2021, 12, 7167–7176 CrossRef CAS PubMed.
- D. E. Smith, H. P. Babcock and S. Chu, Single-Polymer Dynamics in Steady Shear Flow, Science, 1999, 283, 1724–1727 CrossRef CAS PubMed.
- P. LeDuc, C. Haber, G. Bao and D. Wirtz, Dynamics of individual flexible polymers in a shear flow, Nature, 1999, 399, 564–566 CrossRef CAS PubMed.
- S. Gerashchenko and V. Steinberg, Statistics of Tumbling of a Single Polymer Molecule in Shear Flow, Phys. Rev. Lett., 2006, 96, 038304 CrossRef PubMed.
- M. Ripoll, R. G. Winkler and G. Gompper, Star polymers in shear flow, Phys. Rev. Lett., 2006, 96, 188302 CrossRef CAS PubMed.
- C. C. Huang, R. G. Winkler, G. Sutmann and G. Gompper, Semidilute polymer solutions at equilibrium and under shear flow, Macromolecules, 2010, 43, 10107–10116 CrossRef CAS.
- G. Foffi, F. Sciortino, P. Tartaglia, E. Zaccarelli, F. L. Verso and L. Reatto,
et al., Structural arrest in dense star polymer solutions, Phys. Rev. Lett., 2003, 90, 238301 CrossRef CAS PubMed.
- C. Mayer, E. Zaccarelli, E. Stiakakis, C. N. Likos, F. Sciortino and A. Munam,
et al., Asymmetric caging in soft colloidal mixtures, Nat. Mater., 2008, 7, 780–784 CrossRef CAS PubMed.
- J. Mattson, H. M. Wyss, A. Fernandez-Nieves, K. Miyazaki, Z. Hu and D. R. Reichman,
et al., Soft colloids make strong glasses, Nature, 2009, 462, 83 CrossRef PubMed.
- A. M. Philippe, D. Trzzzolillo, J. Galvan-Myoshi, P. Dieudonné-George, V. Trappe and L. Berthier,
et al., Glass transition of soft colloids, Phys. Rev. E, 2018, 97, 040601(R) CrossRef PubMed.
- S. Gupta, J. Stellbrink, E. Zaccarelli, C. N. Likos, M. Camargo and P. Holmqvist,
et al., Validity of the Stokes-Einstein relation in soft colloids up to the glass transition, Phys. Rev. Lett., 2015, 115, 128302 CrossRef PubMed.
- C. N. Likos, A. Lang, M. Watzlawek and H. Löwen, Criterion for determining clustering versus reentrant melting behavior for bounded interaction potentials, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2001, 63, 031206 CrossRef CAS PubMed.
- A. A. Louis, P. G. Bolhuis, J. P. Hansen and E. J. Meijer, Can polymer coils be modeled as “soft colloids”?, Phys. Rev. Lett., 2000, 85, 2522–2525 CrossRef CAS PubMed.
- I. O. Götze, H. M. Harreis and C. N. Likos, Tunable effective interactions between dendritic macromolecules, J. Chem. Phys., 2004, 120, 7761–7771 CrossRef PubMed.
- A. Lang, C. N. Likos, M. Watzlawek and H. Löwen, Fluid and solid phases of the Gaussian core model, J. Phys.: Condens. Matter, 2000, 12, 5087–5108 CrossRef CAS.
- S. Prestipino, F. Saija and P. V. Giaquinta, Phase diagram of the Gaussian-core model, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2005, 71, 050102(R) CrossRef PubMed.
- A. Ikeda and K. Miyazaki, Glass Transition of the Monodisperse Gaussian Core Model, Phys. Rev. Lett., 2011, 106, 015701 CrossRef PubMed.
- A. Ikeda and K. Miyazaki, Slow dynamics of the high density Gaussian core model, J. Chem. Phys., 2011, 135, 054901 CrossRef PubMed.
- D. Coslovich, A. Ikeda and K. Miyazaki, Mean-field dynamic criticality and geometric transition in the Gaussian core model, Phys. Rev. E, 2016, 93, 042602 CrossRef PubMed.
- T. S. Grigera, A. Cavagna, I. Giardina and G. Parisi, Geometric approach to the dynamic glass transition, Phys. Rev. Lett., 2002, 88, 055502 CrossRef PubMed.
- J. M. Bomont, C. N. Likos and J. P. Hansen, Glass quantization of the Gaussian core model, Phys. Rev. E, 2022, 105, 024607 CrossRef CAS PubMed.
- D. Coslovich, M. Bernabei and A. J. Moreno, Cluster glasses of ultrasoft particles, J. Chem. Phys., 2012, 137, 184904 CrossRef PubMed.
- D. Coslovich and A. Ikeda, Cluster and reentrant anomalies of nearly Gaussian core particles, Soft Matter, 2013, 9, 6786–6795 RSC.
- R. Miyazaki, T. Kawasaki and K. Miyazaki, Cluster glass transition of ultrasoft-potential fluids at high density, Phys. Rev. Lett., 2016, 117, 165701 CrossRef PubMed.
- R. Miyazaki, T. Kawasaki and K. Miyazaki, Slow dynamics coupled with cluster formation in ultrasoft-potential glasses, J. Chem. Phys., 2019, 150, 074503 CrossRef PubMed.
- L. Berthier and G. Biroli, Theoretical perspective on the glass transition and amorphous materials, Rev. Mod. Phys., 2011, 83, 587 CrossRef CAS.
- P. Chaudhuri, L. Berthier and W. Kob, Universal nature of particle displacements close to glass and jamming transitions, Phys. Rev. Lett., 2007, 99, 060604 CrossRef PubMed.
- R. Pastore, M. P. Ciamarra, G. Pesce and A. Sasso, Connecting short and long time dynamics in hard-sphere-like colloidal glasses, Soft Matter, 2015, 11, 622–626 RSC.
- K. Kim, Effects of pinned particles on the structural relaxation of supercooled liquids, Europhys. Lett., 2003, 61, 790 CrossRef CAS.
- L. Berthier and W. Kob, Static point-to-set correlations in glass-forming liquids, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2012, 85, 011102 CrossRef PubMed.
- S. Karmakar and G. Parisi, Random pinning glass model, Proc. Natl. Acad. Sci. U. S. A., 2013, 110, 2752–2757 CrossRef CAS PubMed.
- S. Chakrabarty, R. Das, S. Karmakar and C. Dasgupta, Understanding the dynamics of glass-forming liquids with random pinning within the random first order transition theory, J. Chem. Phys., 2016, 145, 034507 CrossRef PubMed.
- W. Kob and L. Berthier, Probing a Liquid to Glass Transition in Equilibrium, Phys. Rev. Lett., 2013, 110, 245702 CrossRef PubMed.
- P. Mausbach and H. O. May, Static and dynamic anomaies in the Gaussian core model liquid, Fluid Phase Equilib., 2006, 249, 17–23 CrossRef CAS.
- H. H. Wensink, H. Löwen, M. Rex, C. N. Likos and S. van Teeffelen, Long-time self-diffusion for Brownian Gaussian-core particles, Comput. Phys. Commun., 2008, 179, 77 CrossRef CAS.
- F. H. Stillinger, Duality relations for the Gaussian core model, Phys. Rev. B: Condens. Matter Mater. Phys., 1979, 20, 299–302 CrossRef CAS.
- P. Charbonneau and G. Tarjus, Decorrelation of the static and dynamic length scales in hard-sphere glass formers, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2013, 87, 042305 CrossRef PubMed.
- M. Z. Slimani, P. Bacova, M. Bernabei, A. Narros, C. N. Likos and A. J. Moreno, Cluster glasses of semiflexible ring polymers, ACS Macro Lett., 2014, 3, 611 CrossRef CAS PubMed.
- L. O. Hedges, L. Maibaum, D. Chandler and J. P. Garrahan, Decoupling of exchange and persistence times in atomistic models of glass formers, J. Chem. Phys., 2007, 127, 211101 CrossRef PubMed.
- R. Pastore, A. Coniglio, A. de Candia, A. Fierro and M. Pica Ciamarra, Cage-jump motion reveals universal dynamics and non-universal structural features in glass forming liquids, J. Stat. Mech., 2016, 054050 CrossRef.
- M. Pica Ciamarra, R. Pastore and A. Coniglio, Particle jumps in structural glasses, Soft Matter, 2016, 12, 358 RSC.
- V. Sposini, A. V. Chechkin, I. M. Sokolov and S. Roldán-Vargas, Detecting temporal correlations in hopping dynamics in Lennard-Jones liquids, J. Phys. A: Math. Theor., 2022, 55, 324003 CrossRef.
- Y. Lanoiselée and D. S. Grebenkov, Unraveling intermittent features in single-particle trajectories by a local convex hull method, Phys. Rev. E, 2017, 96, 022144 CrossRef PubMed.
- P. Mausbach and H. O. May, Transport Anomalies in the Gaussian Core Model Fluid, Z. Phys. Chem., 2009, 223, 1035–1046 CrossRef CAS.
- K. Dawson, G. Foffi, M. Fuchs, W. Götze, F. Sciortino and M. Sperl,
et al., Higher-order glass-transition singularities in colloidal systems with attractive interactions, Phys. Rev. E: Stat., Nonlinear, Soft Matter Phys., 2000, 63, 011401 CrossRef PubMed.
- A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown and P. S. Crozier,
et al., LAMMPS-a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales, Comput. Phys. Commun., 2022, 271, 108171 CrossRef CAS.
-
D. Frenkel and B. Smit, Molecular Simulation: From Algorithms to Applications, Academic Press, New York, 2000 Search PubMed.
-
J. P. Hansen and I. R. McDonald, Theory of Simple Liquids with Applications to Soft Matter, Academic Press, Amsterdam, 4th edn, 2013 Search PubMed.
|
This journal is © The Royal Society of Chemistry 2023 |
Click here to see how this site uses Cookies. View our privacy policy here.