Cluster and reentrant anomalies of nearly Gaussian core particles

Cluster and reentrant anomalies of nearly Gaussian core particles

Daniele Coslovich Email: Université Montpellier 2, Laboratoire Charles Coulomb UMR 5221, Montpellier, France CNRS, Laboratoire Charles Coulomb UMR 5221, Montpellier, France    Atsushi Ikeda Université Montpellier 2, Laboratoire Charles Coulomb UMR 5221, Montpellier, France CNRS, Laboratoire Charles Coulomb UMR 5221, Montpellier, France
July 23, 2019

We study through integral equation theory and numerical simulations the structure and dynamics of fluids composed of ultrasoft, nearly Gaussian particles. Namely, we explore the fluid phase diagram of a model in which particles interact via the generalized exponential potential , with a softness exponent slightly larger than 2. In addition to the well-known anomaly associated to reentrant melting, the structure and dynamics of the fluid display two additional anomalies, which are visible in the isothermal variation of the structure factor and diffusivity. These features are correlated to the appearance of dimers in the fluid phase and to the subsequent modification of the cluster structure upon compression. We corroborate these results through an analysis of the local minima of the potential energy surface, in which clusters appear as much tighter conglomerates of particles. We find that reentrant melting and clustering coexist for softness exponents ranging from up to values relevant for the description of amphiphilic dendrimers, i.e., . .

61.43.Fs, 61.20.Lc, 64.70.Pf, 61.20.Ja
thanks: corresponding author

I Introduction

Hard colloidal suspensions have been extensively studied as the simplest model of particulate systems. In contrast, several colloidal suspensions are composed of soft, macroscopic particles, and display a much more complex behavior Likos_2006 (); Vlassopoulos_Cloitre_2012 (). Star polymers Likos1998 (); zaccarelli_tailoring_2005 (); stiakakis_slow_2010 () and microgels Saunders_Vincent_1999 (); Senff_Richtering_1999 (); Lietor-Santos_Sierra-Martin_Vavrin_Hu_Gasser_Fernandez-Nieves_2009 () are notable examples of purely repulsive, soft colloids, in which softness can be tuned by either altering the macromolecular architecture or by exploiting the sensitivity to temperature and solvent properties. As a consequence these macromolecular aggregates can be compressed to attain very large packing fractions—significantly above random close packing mattsson_soft_2009 (); Koumakis_Petekidis_2012 ().

In theory and simulation, soft colloids can be studied by modeling them as point particles interacting through a soft repulsive potential . Paradigmatic examples are Hertzian particles Pamies_Cacciuto_Frenkel_2009 (), whose repulsion originates from the elastic deformation of spheres at contact, and the Gaussian core model (GCM) Stillinger1976 (); Stillinger1979 (). The thermodynamic and dynamic behaviors of fluids composed of such soft particles differ significantly from the usual hard sphere model and display an anomalous, reentrant density dependence Stillinger1976 (); Stillinger1979 (); louis_mean-field_2000 (); lang_fluid_2000 (); prestipino_phase_2005 (); Mausbach2006 (); mayer_asymmetric_2008 (); krekelberg_anomalous_2009 (); krekelberg_generalized_2009 (); Pamies_Cacciuto_Frenkel_2009 (); Jacquin_Berthier_2010 (); berthier_increasing_2010 (); ikeda_glass_2011 (); Ikeda_therm_2011 (); ikeda_slow_2011 (). Upon increasing density from the dilute regime, the particles first pack more densely and the structural order increases—as in hard spheres. Upon further compression, particles tend to decorrelate spatially, since the energetic gain of avoiding overlaps is overwhelmed by the increase of available configurations (entropic effect) krekelberg_anomalous_2009 (); Jacquin_Berthier_2010 (). This leads to looser spatial correlations and enhanced diffusivity. From a physical point of view, the anomalous behavior can thus be understood in terms of the competition between entropic and energetic effects. The anomalous behavior of the fluid has a direct counterpart in the thermodynamic phase diagram, which displays reentrant melting upon increasing density prestipino_phase_2005 (); Pamies_Cacciuto_Frenkel_2009 (); Ikeda_therm_2011 (). Very recently, such anomalies have indeed been observed in experiments on microgels, which can be modeled rather well using Hertzian particles Paloli_Mohanty_Crassous_Zaccarelli_Schurtenberger_2013 ().

On the other hand, a different type of anomaly can occur if the pair potential is bounded at the origin, i.e., is finite, but has a steep hill at some length scale . An ideal example is the penetrable sphere model (PSM) in which is a finite constant at and zero otherwise Marquest_1989 (). At sufficiently high density, such ultrasoft particles display cluster formation. In the high-density, fluid phase, the radial distribution function shows a clear peak at , meaning that the particles can interpenetrate and form relatively stable clusters Likos_PSM_1998 (). Depending on the details of the interaction potential, formation of cluster crystals Likos_PSM_1998 (); mladek_formation_2006 (); mladek_clustering_2007 (); Zhang_Charbonneau_Mladek_2010 (); Zhang_PSM_2012 () and cluster glasses coslovich_cluster_2012 () have also been observed at high density. In such cluster phases, the collective structure is essentially frozen, but individual particles may diffuse by hopping on the cluster sites moreno_diffusion_2007 (); coslovich_hopping_2011 (). Through a mean-field analysis, Likos et al. Likos_Lang_Watzlawek_Lowen_2001 (); likos_why_2007 () have proposed that clustering occurs whenever the Fourier transform of potential, , has some negative components ( potentials).

Both the reentrant anomalies and clustering can be studied within a coherent model system, called the generalized exponential model (GEM)


where is a softness exponent mladek_formation_2006 (). Clearly, the GCM and the PSM are recovered when and , respectively. The GEM has been extensively studied for and , but very little is known about the behavior for other values of . According to the mean-field argument of Likos et al. Likos_Lang_Watzlawek_Lowen_2001 (); likos_why_2007 (), clustering occurs for . Thus, the model allows to interpolate continuously between reentrant anomalies and cluster formation, with an interesting crossover expected around . Moreover, values of around 3 are of interest, as they may be used to describe the effective interactions between macromolecules with low fractal dimensions, such as amphiphilic dendrimers mladek_computer_2008 (). Although these effective interactions are derived in the infinite dilute regime, a recent study indicates that their phenomenology may be observed in fully atomistic simulation of dendrimers at finite concentrations Lenz_Cluster_2012 ()

In this paper we present a theoretical and numerical investigation of the GEM- with varying softness . We show that the structural and dynamic properties display both reentrant and cluster anomalies when is larger than but close to 2. First, we conduct a preliminary investigation of the fluid phase diagram by means of integral equation theory for various values of . We focus then on the behavior of the model with . We find that the fluid structure and dynamics display a sequence of anomalies, which manifests themselves as multiple minima and maxima in the isothermal properties. We find that the origin of the anomalies observed at high density are due to formation of clusters in the fluid. To provide a sharp characterization of cluster formation we analyze the potential energy surface of the model, following numerical techniques well-established in the context of supercooled liquids sastry__1998 (). An analysis of the local minima sampled above the melting temperature allows to obtain a coherent picture of the structural and dynamic anomalies visible at high density.

Ii Methods

We consider a model of point particles interacting via the generalized exponential potential of order , Eq. (1). As usual, , , and provide the units of length, energy, and time, respectively. We use integral equation theory to predict the static properties over a broad range of state parameters and values. The Ornstein-Zernike equation with the hypernetted chain closure (HNC) relation hansen_theory_2006-1 () was solved by an iterative algorithm. We also used the mean spherical approximation (MSA) closure relation louis_mean-field_2000 (); likos_why_2007 (), which can be solved analytically. Molecular dynamics simulations were carried for systems composed of particles enclosed in a cubic box with periodic boundary conditions. Fluid samples were prepared starting from a random configuration using a Berendsen thermostat to adjust the temperature during the equilibration phase. At low temperature, additional runs were carried out using a slow annealing to make sure that static and dynamic properties were reproducible and were not affected by partial crystallization or phase separation. State points for which at least one sample crystallized were not considered in our analysis, which only focused on the fluid portion of the phase diagram. Production runs were performed in the microcanonical ensemble using the velocity-Verlet algorithm. The integration time step was 0.05 in reduced units for both equilibration and production runs. The properties of the GCM are known to be sensitive to the cut off radius of the potential ikeda_slow_2011 (); ingebrigtsen_what_2012 (). Therefore, a relatively long range cut off was applied (). In addition, we employed a cubic spline interpolation grigera_geometric_2002-1 () between to to improve energy conservation during the MD simulations and to avoid numerical difficulties in the analysis of the potential energy surface. The local minima (or inherent structures stillinger__1982 ()) of the potential energy surface were located by minimizing the total potential energy starting from independent, instantaneous fluid configurations. Thanks to the large system size, a relatively small number of independent configurations (20–40) was sufficient to obtain good quality data. Some larger data sets (up to 300 configurations) were considered at specific state points. The LBFGS algorithm liu__1989 () was used to locate the local minima. We found that convergence becomes particularly slow close to the clustering crossover (see Sec. III). Depending on density and temperature, from 1000 to 50000 iterations were necessary to reduce the norm of the total force per particle to values smaller than . Such a tolerance criterion was sufficient to reduce to zero the number of unstable modes of the dynamical matrix.

Iii Results

The aim of this work is to investigate the interplay between reentrant anomalies and clustering in nearly Gaussian particles, i.e. for values of slightly larger than 2. To understand how the phase diagram is modified as the interactions deviate from the pure Gaussian potential, we start with a preliminary exploration of fluid structure of the model by varying systematically the exponent . For this purpose, we employ the integral equation theory using the HNC closure relation, which provides a good description of the structure of soft particles at sufficiently high density likos_effective_2001 (); mladek_clustering_2007 ().

A simple way to highlight the anomalous features of the structure in soft core particles is to identify the locus of state points in the plane for which the height of the first peak of the structure factor is constant. We found that HNC suffers from convergence problems at intermediate densities, in particular close to the clustering crossover (see below), when . Therefore, we restrict ourselves to a relatively small value of . For the models at hand, we use the condition , which is smaller than typical values observed at melting in simple liquids ( according to the Hansen-Verlet criterion Hansen_Verlet_1969 ()). We emphasize that the aim of this preliminary analysis is to highlight the presence of a structural anomaly in the diagram and not to track closely the fluid-solid phase boundary. Such an operation would require adjusting the value of depending on density, since the first peak of reaches large values (i.e., ) along the cluster crystallization line likos_why_2007 (); mladek_clustering_2007 ().

The iso- curves obtained within the HNC approximation are shown in Fig. 1 for several values of , along with the fluid-solid phase boundary of the GCM determined in Refs. prestipino_phase_2005 () and Ikeda_therm_2011 (). At low density, HNC predicts the iso- line to depend very little on and to follow the low density branch of the GCM phase boundary. In this portion of the phase diagram, the HNC predictions are semi-quantitative at best, since we expect the Hansen-Verlet criterion to work well along the low density phase boundary. Thus, larger values of would be expected. For all the values of shown in the figure, the iso- lines show a first decrease upon increasing density, following the non-monotonic, anomalous behavior of the melting line of the GCM. Upon increasing the density even further, an additional branch appears: for values of strictly larger than 2, the iso- lines follow a common envelope until they reach a crossover density and start bending upwards. The crossover density diverges when , as it can be inferred by tracking the location of the minimum of as a function of (not shown). Such a behavior is a striking confirmation of the mean-field argument of Likos et al. Likos_Lang_Watzlawek_Lowen_2001 (). In fact, as it will be clear from the following, the minimum of the iso- lines marks the onset of clustering in the fluid phase. We find that the two coexisting anomalies are visible up to values of around 3, at which the minimum transforms into a small inflection. Above the iso- lines predicted by HNC display a purely monotonic increase with density.

Figure 1: Integral equation theory predictions for the iso- lines in the - plane for various values of . The lines connect state points for which the height, , of the first peak of the structure factor is equal to 2.3. Full and dashed lines correspond to HNC and MSA approximations, respectively. Open squares indicate the fluid-solid phase boundary of the GCM () as determined in Refs. prestipino_phase_2005 () and  Ikeda_therm_2011 (). For clarity, the iso- lines predicted by MSA are not shown below the fluid-solid phase boundary.

Results obtained using the MSA are also included in the figure. In this case, the iso- lines can be obtained analytically and are given by straight lines in the - plane

where is the value of the minimum of the Fourier transform of the potential. We see that the MSA predictions follow closely those obtained within HNC at high density. However, the crossover between the low and intermediate density regimes (reentrant anomaly) is not captured within the mean-field description of MSA, as the latter is only reliable at high temperatures and densities.

For values of the softness exponents between and approximately 3, HNC thus predicts the existence of three regimes in the fluid phase diagram separated by two anomalies: (i) the “reentrant” anomaly well known from studies of the GCM and Hertzian spheres, corresponding to a maximum in the iso- lines, and (ii) a “cluster” anomaly, observed as the particles start to explore the ultrasoft portion of the potential, which appears as a minimum in the iso- lines. To test this anomalous scenario in-depth, we focus in the following on the case , for which we collected the most extensive set of simulation data. Simulations for different values of will be briefly discussed at the end of this section.

Figure 2: Fluid phase diagram of the GEM-2.2 in the - plane. State points at which are indicated as full squared and are connected by full lines (iso- lines). State points at which are drawn as filled or open circles according to whether they are characterized by clustering or transient overlaps, respectively, as explained in the text. Open squares indicate the fluid-solid phase boundary of the GCM () as determined in Refs. prestipino_phase_2005 () and  Ikeda_therm_2011 (). The light grey and dark grey shaded areas indicate portions of the phase diagram displaying finite fractions of dimers and trimers, respectively, in the local minima of the potential energy surface. Inset: magnification of the behavior of the iso- line at large densities for (squares) and (triangles). Also included are the corresponding results from HNC (solid lines).

In Fig. 2 we combine our molecular dynamics results and the analysis of the potential energy surface to characterize the fluid phase diagram the GEM-2.2. We start by discussing the shape of the iso- line, which we determined using the value (Hansen-Verlet criterion). At low density, the line first follows the trend of the GCM fluid-solid phase boundary. Upon further compression, it first passes through a maximum located around and then a minimum around . Thus both the reentrant and cluster anomalies predicted by HNC are visible in the simulation data. Upon closer inspection, however, it is clear that the agreement between HNC and the simulations in this portion of the phase diagram is only qualitative. The minimum separating the soft-core and the cluster regimes is in fact deeper than predicted by HNC. These discrepancies are further highlighted in Figs. 3(a) and (c), where we compare the HNC predictions for the structure factor and the radial distribution function around the crossover. HNC predicts clustering to occur at smaller densities (or higher temperatures) than actually observed in the simulations, as seen from the differences in around . Only at larger densities and temperatures [see Figs. 3(b) and (d)], the HNC predictions agree quantitatively with the simulations, as expected from previous work mladek_clustering_2007 (). This is further illustrated by the iso- lines shown in the inset of Fig. 2, where we observe good agreement between HNC and simulations for , but not for . We conclude that better closure relations are required for a quantitative description of the fluid phase diagram close to the fluid-solid phase boundary, especially at low and intermediate densities.

Figure 3: Comparison of the structure factors [panels (a) and (b)] and radial distribution functions [panel (c) and (d)] obtained using the HNC closure (solid lines) and MD simulations (crosses). Data are shown for in (a) and (c), and for in (b) and (d).

To illustrate the connection between the minimum of the iso- line and clustering, we evaluate the radial distribution function and analyze its isothermal variation. In the following, we will denote with the position of the peak of corresponding to the first neighbor (FN) shell. In Fig. 4(a) we show the along the isotherm . As the density is increased from , the first peak initially shifts markedly towards smaller distances, as expected for soft-core particles. Upon crossing , a broad peak develops at and grows by further increasing the density. The presence of a well-defined minimum between and , indicates the formation relatively stable clusters. The FN peak position also shows an interesting behavior: it first decreases with increasing density up to , then it starts to increase until the density reaches approximately . This point will be further discussed below.

In contrast to cluster formation, transient overlaps can always occur at any density in fluids of ultrasoft particles as a consequence of thermal fluctuations. Indeed, at high temperatures, coreless particles often superpose during collisions, which implies a finite value of . How to distinguish between clustering and such transient overlaps? We propose to classify clustering states as those for which the second derivative . In practice, an equivalent and numerically more stable criterion is to check whether the displays a minimum between 0 and . States for which is finite but are characterized instead by transient overlaps.

Figure 4: Density variation of the radial distribution function and cluster properties of the GEM-2.2 evaluated for instantaneous fluid configurations: (a) along the isotherm . For clarity, each line has been vertically shifted from the previous one by a constant offset, 1.25. For , the horizontal dotted lines indicate the ordinates at which would equal zero. Filled circles indicate the first neighbor shell peaks. (b) Fraction of -mers for as a function of for (monomers, triangles), 2 (dimers, circles), and 3 (trimers, squares). (c) and (d) Monomer-monomer (dotted line), monomer-dimer (dashed line) and dimer-dimer (solid line) radial distribution functions evaluated at (c) , and (d) , . Note that the distribution functions are zero by construction for .

With all this information at hand we return to the fluid phase diagram Fig. 2 and mark the boundary between normal and cluster states. To this end, we identify the state points in the - plane at which reaches constant, “small” value. The reference value is small enough to detect the early onset of clusters and overlaps, still sufficiently large to minimize statistical uncertainties. Numerically, we evaluate from the probability of finding particles at distances less than 0.03. Following the classification introduced above, we use open or filled circles in Fig. 2 for states characterized by clusters or transient overlaps, respectively. In the low and intermediate density regime, the onset of transient overlaps is located at temperatures well above the melting line. Upon increasing the density, the iso- line bends downwards and starts to run vertically in the plane. In this portion of the phase diagram, the sign of along iso- line becomes negative, signaling the transition to cluster states. By visual inspection, it is clear that the boundary between normal and cluster states corresponds well to the minimum of the iso- line, which motivates the name “cluster anomaly”.

To understand the physical origin of the cluster anomaly, we analyze in more detail the changes in the microscopic structure of the clusters. To identify clusters, we use the same geometrical definition as in previous work coslovich_hopping_2011 (); coslovich_cluster_2012 (): a particle is considered to be part of a cluster if its distance to at least another particle of that cluster is smaller than a fixed value . Typically, is chosen around the first minimum of , which is located in the range 0.6–0.7 over the pertinent range of density and temperatures. For simplicity, in our analysis we used a state independent cut off .

The cluster population , where is the number of -mers in a configuration, is shown in Fig. 4(b) as a function of density along the isotherm . We find that starting from the fluid is no longer composed only of isolated particles (monomers) but possesses a growing fraction of dimers. In the range of density between 0.78 and 1.1, we observe coexistence between monomers and dimers, the fraction of higher order -mers being zero within our numerical accuracy. We remark that this range varies slightly as a function of the chosen temperature. To illustrate the real space structure of coexisting monomers and dimers, we show in Fig 4(c) and (d) the monomer-monomer (MM), monomer-dimer (MD) and dimer-dimer (DD) radial distribution functions evaluated at and , respectively. These correlation functions are calculated from

where and are either or and the term in the double sum is ignored when . For dimers, the distance is evaluated from their center of mass. We find that monomers are more loosely correlated than dimers, as it can be evinced from the nearly flat shape of . This probably reflects the fact the effective interactions between dimers are harsher than the bare potential between individual particles. Dimers, on the other hand, tend to occupy a larger volume and are thus separated from their neighbors by larger distances. This suggests that the slight increase of above is due to the larger FN shell of the growing dimers population.

Upon further compression above , a small fraction of trimers appears and grows by increasing density at the expenses of and [see Fig. 4(b)]. Keeping fixed, we cannot extend the analysis of the fluid regime beyond , due to incipient, partial cluster crystallization. We thus conclude that, above , clusters are characterized by low population numbers, i.e. 2 or 3, and coexist in a nearly structureless continuum formed by isolated particles. This scenario is most likely connected to the underlying crystalline phase diagram, which may possess a cascade of transitions involving normal fluids (i.e., isolated particles) and cluster crystals with integer population numbers. Such a cascade of complex crystal structures is known to occur in the PSM Likos_PSM_1998 (); Zhang_PSM_2012 () and in the GEM-4 Zhang_Charbonneau_Mladek_2010 (); Neuhaus_Likos_2011 ().

Figure 5: Density variation of the radial distribution function and cluster properties of the GEM-2.2 evaluated for local minima of the potential energy surface. Lines and symbols as in Fig. 4. Note that in this case both the radial distribution function and the cluster populations are independent of over the studied range of state parameters. The actual data shown in the figure were obtained from minimizations at . In (a), solid lines are colored according to the different regimes of cluster populations, also indicated by the shaded areas in (b). Vertical arrows highlight the splitting of the FN maxima at and and indicate the location of the secondary peaks. In (c) and (d) the partial radial distribution functions are shown for , and , , respectively.

Due to the thermal motion of particles, however, the definition of a cluster in the fluid states is subject to some ambiguity. For instance, the drop of at for (and to a less extent for ) indicates that there is no clear separation between the monomers and the FN shells of the clusters. To obtain a sharper picture, we apply our cluster analysis to the local minima of the underlying potential energy surface (PES). We perform a statistical sampling of the PES properties by quenching several independent samples (instantaneous configurations) at a given and , as described in the Sec. II. In the following, we will focus on the isotherm and study how the underlying landscape changes upon varying the density. Although the PES is fully determined once the number of particles and the volume are specified, the properties of the basins visited by the system may vary below some onset temperature sastry__1998 (). However, this is not observed for the GEM-2.2, at least down to temperatures close to the fluid-solid phase boundary. The situation resembles the scenario observed in simple liquids, where the properties of local minima visited above melting are indeed temperature independent sastry__1998 ().

The main advantage over the analysis based on instantaneous configurations is that, in the local minima, clusters appear as much tighter conglomerates of particles, sitting practically on top of each other. This is demonstrated by the radial distribution function evaluated using the minimized configurations, see Fig. 5(a). The peak at zero distance, which emerges at the clustering crossover, is now extremely narrow and is separated from the FN shell by a clear gap where is zero within numerical precision. This, in turn, makes the definition of clusters straightforward and unambiguous. In practice, we used as a cut off distance for the cluster analysis on local minima.

Compared to instantaneous configurations, local minima present a more subtle density dependence of the radial distribution function. The latter is best understood when discussed along with that of the cluster population , which we show in Fig. 5(b). An analysis of the cluster populations reveals a steep rise of dimers at , making the definition of the clustering crossover very sharp. The fraction of trimers becomes non zero only above , again displaying a sharper growth than in instantaneous configurations. Upon increasing the density even further (), monomers eventually disappear and we observe coexistence between dimers and trimers. These various regimes are indicated in Fig. 2 by shaded areas, to highlight their connection to the observed anomalies.

We can now more easily understand the evolution of with density. The FN peak of shows a broadening around and then splits upon further increasing the density. The two peaks arise from the coexistence of monomers and dimers, which are characterized by two distinct typical inter-cluster separations. This is confirmed by an analysis of the monomer-monomer and dimer-dimer contributions to , which are calculated as for instantaneous configurations and are shown in Fig. 5(c) and (d) at and , respectively. In particular, the splitting of the FN shell for [see panel (a)] is due to growing dimer-dimer correlations. Indeed, the position of the secondary maximum at matches the location of the first peak in . We also find that the MM peak at small distances soon looses its structure above , whereas the DD peak shifts progressively to smaller distances. This, together with the density evolution of above , suggests that dimers behave as larger, effective soft-core particles. Above , the FN peak of shows again a broadening and a subsequent splitting, which arises now from the coexistence of subpopulations of dimers and trimers in the local minima.

The partial radial distribution functions for [see Fig. 5(d)] qualitatively confirm the results obtained from instantaneous configurations in the monomer-dimer coexistence region. In particular, we see that isolated particles, which amount to of the -mers at this density, are spatially uncorrelated, as indicated by the unstructured profile of . Very close to the cluster crossover [see Fig. 5(c)], we observe a strong suppression of and and a concomitant enhancement of monomer-monomer correlations. This indicates that monomers and dimers have a tendency to microsegregate and suggests an underlying instability with respect to phase separation into coexisting fluid and cluster phases. Visual inspection of the real space structure of the local minima at confirms these observations. We found no such tendency to microsegregation in the instantaneous configurations. These results resemble the scenario observed in certain phase separating binary Lennard-Jones mixtures sarkar_inherent_2011 (), for which the local minima display clear signs of demixing—even for states at which the equilibrium fluid is essentially homogeneous.

Figure 6: Structural and dynamic anomalies of the GEM-2.2 as a function of density along selected isotherms: (squares), 0.014 (circles), 0.012 (triangles), 0.011 (reversed triangles), and 0.01 (diamonds). (a) Height, , of the first peak of the structure factor. (b) Inverse of two-body excess entropy . (c) Diffusion constant . Arrows at and mark the appearance of dimers and trimers in the PES, respectively.

We now analyze the connection between the structural anomalies and the dynamics. It is well established krekelberg_anomalous_2009 () that the reentrant anomaly in soft-core models is accompanied by a dynamic anomaly in the single particle diffusivity: beyond the reentrant crossover, the dynamics accelerates upon increasing density due to the progressive the loss of structural correlations. We now show that also the cluster anomalies also have their dynamic counterparts. In the three panels of Fig. 6 we compare the isothermal variation of the peak height of the structure factor, the inverse of the two-body contribution to the excess entropy

and the diffusion constant , which has been obtained from the usual Einstein relation

where is the mean square displacement. The two-body excess entropy is a simple measure of pair translation order in the system. This quantity has been used to characterize structural anomalies in fluids of ultrasoft particles krekelberg_anomalous_2009 (); krekelberg_generalized_2009 ().

Figure 7: Structural and dynamic anomalies of the GEM-3 as a function of density along selected isotherms: (squares), 0.05 (circles), 0.04 (triangles), 0.03 (reversed triangles), and 0.02 (diamonds). (a) , of the first peak of the structure factor. (b) Inverse of the two-body excess entropy . (c) Diffusion constant .

The two isotherms at high temperature reproduce the reentrant anomalies in both structure and dynamics, as expected. On the contrary the two low temperature isotherms necessarily hit, at small , on the fluid-bcc phase boundary, and therefore only the cluster anomalies are visible. We find that follows closely the density variation of the . A close inspection of panels (a) and (b) of Fig. 6 reveals only a minor discrepancy between the location of the minima in and , with the latter providing a slightly larger () value of the density at the minimum. A comparison of the structural metrics and shows that the clustering crossover at is clearly accompanied by a decrease in diffusivity. Physically, this anomaly can be easily understood, since transient clusters temporarily behave as effective, heavier particles, thus slowing down the dynamics of the fluid. Despite the clear connection between structure and dynamics illustrated by Fig. 6, an attempt to collapse the dynamic data as a function of the two-body excess entropy , as in previous works krekelberg_anomalous_2009 (); krekelberg_generalized_2009 (), did not yield satisfactory results on a quantitative level (see also Fomin_Ryzhov_Gribova_2010 () for a discussion on the breakdown of such a “Rosenfeld scaling” rosenfeld_quasi-universal_1999 ()).

Figure 8: Fraction of -mers in local minima of the potential energy surface as a function of for (triangles), 2 (circles), 3 (squares), and 4 (reversed triangles). Local minima were sampled at , .

Around and at sufficiently low temperatures, an additional anomaly is visible as a maximum in (and ) and a minimum in . The relatively large values attained by along the lowest isotherm, , suggest that the system might be slightly supercooled, despite being (meta)stable during the whole length of our simulations. Figure 6(a) and (b) show however that this maximum develops smoothly by decreasing temperature. The inflection of the iso- lines in the diagram [see inset of Fig. 2] is most likely a high temperature manifestation of this additional anomaly. Thus, we are confident that this represents an actual feature of the equilibrium fluid. This anomalous scenario is also observed at smaller values, down to at least 2.01 (not shown). In this case, the two cluster anomalies are visible at much lower temperatures than in the GEM-2.2, as expected from the HNC results.

Qualitatively, the second cluster anomaly in the GEM-2.2 correlates with the appearance of trimers and the disappearance of isolated particles. However, its underlying physical origin remains unclear: for instance, why should the formation of trimers accelerate the dynamics? A possible explanation is that, at these densities, different subpopulations of -mers coexist (), which enhances the structural disorder and thus accelerates the dynamics. It could also be argued that dimers behave as effective, ultrasoft particles and display an additional reentrant anomaly in their own right. We note, in passing, that neither HNC nor MSA can reproduce this anomalous feature. Further work is definitely needed to clarify this issue.

By varying the softness exponent, we found that coexistence of reentrant and cluster anomalies is well visible up to at least , a value that can be used to model the effective interactions between amphiphilic dendrimers mladek_computer_2008 (). This is demonstrated in Fig. 7, where we show the isothermal variation of , and in the GEM-3. In this case, however, only the first cluster anomaly is visible in the range of densities and temperatures covered by our simulations. The absence of the second anomaly is consistent with the fact the fraction of trimers in the fluid is negligible () at the state points shown in Fig. 7.

Finally, we study the evolution of the cluster population in the local minima as . Figure 8 shows the fraction of -mers with at fixed and . As the interaction smoothly transforms into a Gaussian, clusters of high order disappear progressively. At the studied density, only models with are fully dissociated within our statistics. The variation of with thus provides a striking confirmation of the Likos criterion for clustering in fluid states.

Iv Conclusions

In this work we demonstrated the existence of multiple anomalies in the structure and dynamics of fluids composed of nearly Gaussian particles. We focused in particular on particles interacting via the generalized exponential potential, Eq. (1), with an exponent slightly larger than 2. Our most extensive set of data concerns the GEM-2.2, for which we gathered molecular dynamics results and for which we analyzed systematically the potential energy surface.

At low densities, the structure and dynamics of the fluid display a non-monotonic density variation. Such a reentrant anomaly is well-known in soft-core particles krekelberg_anomalous_2009 (); Pamies_Cacciuto_Frenkel_2009 () and manifests itself as a loss of translation order upon compression and a concomitant acceleration of the single-particle dynamics. Nearly Gaussian particles show two additional anomalies by compressing the fluid at constant temperature: the height of the first peak of the structure factor displays a minimum and an additional inflection at higher densities. At sufficiently low temperatures, such an inflection transforms into a small but well-defined minimum. All these anomalies have counterparts in the single particle dynamics, with minima of corresponding to maxima of the diffusion constant , and viceversa. We note that a simple scaling between Newtonian and overdamped Langevin dynamics has been recently reported for ultrasoft fluids Pond_Errington_Truskett_2011 (), indicating that the connection between structural and dynamical anomalies observed here is likely quite general.

At high densities, particles explore the bounded part of the potential. The two additional anomalies occurring in this portion of the fluid phase diagram, can thus be attributed to the appearance of clusters in the fluid. Namely, we found that the minimum in correlates with the sharp growth of dimers in both instantaneous configuration and in the local minima of the PES. The additional maximum ostensibly correlates to the appearance of trimers upon further compression, although its actual physical origin needs to be clarified. Interestingly, the spatial structure of the local minima reveals that coexisting subpopulations of isolated particles and clusters can show a slight segregation. This suggests an underlying instability of fluid with respect to phase separation into coexisting fluid and cluster phases. It would thus be interesting to characterize in more detail the low temperature phase diagram of nearly Gaussian particles, along the lines of recent numerical work on the GEM-4 Zhang_Charbonneau_Mladek_2010 ().

We demonstrated that clustering in the fluid phase disappears in a continuous way as the interaction potential tends to a purely Gaussian shape, i.e., as the softness exponent approaches 2 from above. This is an explicit confirmation of the mean-field argument by Likos et al. Likos_Lang_Watzlawek_Lowen_2001 (), which has never been so far tested explicitly on GEMs with close to 2. What our results clarify, however, is that reentrant and clustering behavior are not mutually exclusive, as is sometimes implicit in the presentation of and models Likos_Lang_Watzlawek_Lowen_2001 (). We found that coexistence of reentrant and cluster anomalies is well visible up to at least . This behavior was anticipated by Zhang et al Zhang_Charbonneau_Mladek_2010 () in their study of the low temperature phase diagram of the GEM-4 model and has been confirmed here more systematically. Whether these coexisting anomalies may be observed in colloidal suspensions of highly branched macromolecules, such as dendrimers, is an open question that awaits an experimental verification.

Another interesting implication of our work concerns the nature of the possible glassy states of ultrasoft particles. It has been shown that the GCM becomes a surprisingly good glass-former at sufficiently high densities (ikeda_glass_2011 (). This has been attributed to a frustration effect induced by the long range of the Gaussian potential. We found that for the cluster anomaly is shifted up to approximately . This, in turn, opens up a broad region in the fluid phase diagram in which the model behaves essentially as the GCM. In such a region, models with values marginally larger than 2 might display peculiar glass-forming properties. In particular, it would be interesting to see whether the proximity to cluster formation destabilizes the glassy behavior of the GCM and to what extent such a behavior will be affected by clustering at even higher densities. Finally, we believe that the analysis of the PES may prove particularly useful in the study of amorphous cluster phases of size dispersed ultrasoft particles coslovich_cluster_2012 (). Work to investigate the properties of the PES explored by cluster glass-forming fluids is currently underway.

We acknowledge the HPC@LR Center of Competence in High-Performance Computing of Languedoc-Roussillon (France) for allocation of CPU time. A. I. acknowledges the financial support from JSPS postdoctral fellowship for research abroad.


  • (1) C. N. Likos, Soft Matter 2, 478 (2006)
  • (2) D. Vlassopoulos and M. Cloitre, Soft Matter 8, 4010 (2012)
  • (3) C. N. Likos, H. Löwen, M. Watzlawek, B. Abbas, O. Jucknischke, J. Allgaier, and D. Richter, Phys. Rev. Lett. 80, 4450 (1998)
  • (4) E. Zaccarelli, C. Mayer, A. Asteriadi, C. N. Likos, F. Sciortino, J. Roovers, H. Iatrou, N. Hadjichristidis, P. Tartaglia, H. L�wen, and D. Vlassopoulos, Phys. Rev. Lett. 95, 268301 (2005)
  • (5) E. Stiakakis, A. Wilk, J. Kohlbrecher, D. Vlassopoulos, and G. Petekidis, Phys. Rev. E 81, 020402 (2010)
  • (6) B. R. Saunders and B. Vincent, Adv. Colloid Interface Sci. 80, 1 (1999)
  • (7) H. Senff and W. Richtering, J. Chem. Phys. 111, 1705 (1999)
  • (8) J.-J. Lietor-Santos, B. Sierra-Martin, R. Vavrin, Z. Hu, U. Gasser, and A. Fernandez-Nieves, Macromolecules 42, 6225 (2009)
  • (9) J. Mattsson, H. M. Wyss, A. Fernandez-Nieves, K. Miyazaki, Z. Hu, D. R. Reichman, and D. A. Weitz, Nature 462, 83 (2009)
  • (10) N. Koumakis, A. Pamvouxoglou, A. S. Poulos, and G. Petekidis, Soft Matter 8, 4271 (2012)
  • (11) J. C. Pamies, A. Cacciuto, and D. Frenkel, J. Chem. Phys. 131, 044514 (2009)
  • (12) F. H. Stillinger, J. Chem. Phys. 65, 3968 (1976)
  • (13) F. H. Stillinger, Phys. Rev. B 20, 299 (1979)
  • (14) A. A. Louis, P. G. Bolhuis, and J. P. Hansen, Phys. Rev. E 62, 7961 (2000)
  • (15) A. Lang, C. N. Likos, M. Watzlawek, and H. Lowen, J. Phys.: Condens. Matter 12, 5087 (2000)
  • (16) S. Prestipino, F. Saija, and P. V. Giaquinta, Phys. Rev. E 71, 050102 (2005)
  • (17) P. Mausbach and H.-O. May, Fluid Phase Equilibria 249, 17 (2006)
  • (18) C. Mayer, E. Zaccarelli, E. Stiakakis, C. N. Likos, F. Sciortino, A. Munam, M. Gauthier, N. Hadjichristidis, H. Iatrou, P. Tartaglia, H. Lowen, and D. Vlassopoulos, Nature Mat. 7, 780 (2008)
  • (19) W. P. Krekelberg, T. Kumar, J. Mittal, J. R. Errington, and T. M. Truskett, Phys. Rev. E 79, 031203 (2009)
  • (20) W. Krekelberg, M. Pond, G. Goel, V. Shen, J. Errington, and T. Truskett, Phys. Rev. E 80, 061205 (2009)
  • (21) H. Jacquin and L. Berthier, Soft Matter 6, 2970 (Jun 2010)
  • (22) L. Berthier, A. J. Moreno, and G. Szamel, Phys. Rev. E 82, 060501 (2010)
  • (23) A. Ikeda and K. Miyazaki, Phys. Rev. Lett. 106, 015701 (2011)
  • (24) A. Ikeda and K. Miyazaki, J. Chem. Phys. 135, 024901 (2011)
  • (25) A. Ikeda and K. Miyazaki, J. Chem. Phys. 135, 054901 (2011)
  • (26) D. Paloli, P. S. Mohanty, J. J. Crassous, E. Zaccarelli, and P. Schurtenberger, Soft Matter 9, 3000 (2013)
  • (27) C. Marquest and T. A. Witten, J. Phys. France 50, 1267 (1989)
  • (28) C. N. Likos, M. Watzlawek, and H. Löwen, Phys. Rev. E 58, 3135 (1998)
  • (29) B. M. Mladek, D. Gottwald, G. Kahl, M. Neumann, and C. N. Likos, Phys. Rev. Lett. 96, 045701 (2006)
  • (30) B. M. Mladek, D. Gottwald, G. Kahl, M. Neumann, and C. N. Likos, J. Phys. Chem. B 111, 12799 (2007)
  • (31) K. Zhang, P. Charbonneau, and B. M. Mladek, Phys. Rev. Lett. 105, 245701 (2010)
  • (32) K. Zhang and P. Charbonneau, J. Chem. Phys. 136, 214106 (2012)
  • (33) D. Coslovich, M. Bernabei, and A. J. Moreno, J. Chem. Phys. 137, 184904 (2012)
  • (34) A. J. Moreno and C. N. Likos, Phys. Rev. Lett. 99, 107801 (2007)
  • (35) D. Coslovich, L. Strauss, and G. Kahl, Soft Matter 7, 2127 (2011)
  • (36) C. N. Likos, A. Lang, M. Watzlawek, and H. Löwen, Phys. Rev. E 63, 031206 (2001)
  • (37) C. N. Likos, B. M. Mladek, D. Gottwald, and G. Kahl, J. Chem. Phys. 126, 224502 (2007)
  • (38) B. M. Mladek, G. Kahl, and C. N. Likos, Phys. Rev. Lett. 100, 028301 (2008)
  • (39) D. A. Lenz, R. Blaak, C. N. Likos, and B. M. Mladek, Phys. Rev. Lett. 109, 228301 (2012)
  • (40) S. Sastry, P. G. Debenedetti, and F. H. Stillinger, Nature 393, 554 (1998)
  • (41) J.-P. Hansen and I. R. McDonald, Theory of Simple Liquids, 3rd ed. (Academic Press, New York, 2006)
  • (42) T. S. Ingebrigtsen, T. B. Schrøder, and J. C. Dyre, Phys. Rev. X 2, 011011 (2012)
  • (43) T. S. Grigera, A. Cavagna, I. Giardina, and G. Parisi, Phys. Rev. Lett. 88, 055502 (2002)
  • (44) F. H. Stillinger and T. A. Weber, Phys. Rev A 25, 983 (1982)
  • (45) D. C. Liu and J. Nocedal, Math. Program. 45, 503 (1989)
  • (46) C. N. Likos, Phys. Rep. 348, 267 (2001)
  • (47) J.-P. Hansen and L. Verlet, Phys. Rev. 184, 151 (1969)
  • (48) T. Neuhaus and C. N. Likos, J. Phys.: Condens. Matter 23, 234112 (2011)
  • (49) S. Sarkar and B. Bagchi, Phys. Rev. E 83, 031506 (2011)
  • (50) Y. D. Fomin, V. N. Ryzhov, and N. V. Gribova, Phys. Rev. E 81, 061201 (2010)
  • (51) Y. Rosenfeld, J. Phys.: Condens. Matter 11, 5415 (1999)
  • (52) M. J. Pond, J. R. Errington, and T. M. Truskett, Soft Matter 7, 9859 (2011)
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description