Colloidal Gelation with Variable Attraction Energy

Colloidal Gelation with Variable Attraction Energy


We present an approximation scheme to the master kinetic equations for aggregation and gelation with thermal breakup in colloidal systems with variable attraction energy. With the cluster fractal dimension as the only phenomenological parameter, rich physical behavior is predicted. The viscosity, the gelation time and the cluster size are predicted in closed form analytically as a function of time, initial volume fraction and attraction energy by combining the reversible clustering kinetics with an approximate hydrodynamic model. The fractal dimension modulates the time evolution of cluster size, lag time and gelation time and of the viscosity. The gelation transition is strongly nonequilibrium and time-dependent in the unstable region of the state diagram of colloids where the association rate is larger than the dissociation rate. Only upon approaching conditions where the initial association and the dissociation rates are comparable for all species (which is a condition for the detailed balance to be satisfied) aggregation can occur with . In this limit, homogeneous nucleation followed by Lifshitz-Slyozov coarsening is recovered. In this limited region of the state diagram the macroscopic gelation process is likely to be driven by large spontaneous fluctuations associated with spinodal decomposition.

I Introduction

Colloidal suspensions gel if there is an attractive interactions of sufficient strength between the particles. This gelation transition has been the focus of intense research during the last decade since it plays an important role in many practical applications as e.g. processing of polymers or food technology. In spite of intensive efforts (1); (2); (3); (4); (5); (6); (7); (8); (9); (10); (11); (12); (13) in the past aimed at clarifying the nature of the gelation transition, the basic mechanism by which a fluid colloidal suspension turns into solid remains unclear. Many numerical studies have been proposed over the last decades which have brought a wealth of phenomenological information about the connection between microscopic attraction and the gelation process (14); (15); (16); (17); (18); (19); (20). However, analytical models are lacking, and therefore it is difficult to elucidate the basic mechanisms and to extract scaling laws in analytical form.

Some time ago, the idea has been proposed (21) that the gelation transition may be interpreted as a ”renormalized” glass transition where the growing colloidal clusters occupy an increasingly larger volume fraction up to the point at which their motions become governed by glassy correlation, the clusters become caged by their neighbors and the system becomes solid by interconnection or random packing of clusters. This scenario is different from what one observes in chemical gels where the bonding is permanent (in contrast with colloidal bonds that can be broken up by thermal energy) and percolation provides an excellent description of chemical gelation (22). With colloidal gels, however, simulations (23) have established that the dynamics is strikingly different from that of chemical gels and colloidal gelation cannot be understood with percolation concepts alone. The concept of colloidal gelation as a cluster-jamming transition has clearly brought progress in the modeling of the static structure-elasticity relation of dense colloidal gels (24). However, it has not been implemented in an analytical model of the gelation transition that can be tested in comparison with experiments. The main problem resides in the difficulty of bridging the macroscopic mechanical response (the viscosity) with the mesoscopic level of the clusters and ultimately with the underlying microscopic association/dissociation kinetics of individual colloidal particles.

Figure 1: (color online). Schematic representation of the thermoresponsive nanoparticles used in the gelation experiments. The dependence of the attraction on temperature (Eq. (31)) is also shown in comparison with the experimental data of Ref. (25). is the highest attraction energy reachable with this system.

Here we present an analytical model of colloidal gelation with variable attraction energy. The model provides a framework which connects the level of the pair-attraction energy with the mesoscopic level of the clusters and finally with the macroscopic mechanical response. Analytical laws can be extracted for the viscosity, the gelation time and the cluster size. These laws provide a theoretical explanation to several observations in the past for which no theoretical description is available.

This investigation was prompted by our finding that well-defined attractive interaction can be induced in suspensions of thermosensitive microgels by raising the temperature above the volume transition of these systems.(25) Figure 1 shows these particles in a schematic fashion: A dense network of the thermosensitive poly(N-isopropylacrylamide) (PNIPAM) is grafted to solid polystyrene core. Immersed in cold water the shell of these particles will swell. Above the volume transition at 32-33 most of the water will be expelled from the network because water becomes a poor solvent under these conditions. As a consequence, the shrunken shells will become mutually attractive and the strength of the ensuing attractive interaction can be adjusted precisely.(25) The advantages of this system for the study of gelation are at hand: Well-defined short-range attraction can be induced without adding an additional component by just raising the temperature. Moreover, the fluid state can be recovered by simply lowering the temperature below 32. The sharp transition observed experimentally can be traced back to the predicted strongly nonlinear dependence of the viscosity on .

This article is structured as follows. First we present the analytical theory of colloidal aggregation and gelation with variable attraction energy. Then we discuss the predictions in terms of scaling laws and a state diagram extracted from the theory. Then we describe the gelation experiments with thermosensitive colloidal particles. Finally, we present the comparison between theory and experiments.

Ii Kinetics of aggregation and gelation with variable attraction energy

ii.1 Assumptions and steps in the derivation

The model is based upon the following steps and assumptions. (i) Any two colloidal particles interact via a rectangular-well attractive potential of width and depth . (ii) The clustering process is described by a master kinetic equation with an effective association rate which accounts for bond dissociation. (iii) The association and dissociation rates between two particles are evaluated from steady-state solutions to the Smoluchowski diffusion equations for the two-body dynamics in a mutual attraction potential. (iv) The so obtained analytical solution of the master kinetic equation is used to obtain analytical expressions of the time-dependent cluster size distribution and of the time-dependent volume fraction occupied by the clusters in the system. (v) Since the attraction is short-range and the hydrodynamics is screened from the clusters interior, clusters are assumed to behave hydrodynamically like hard-spheres and the time evolution viscosity of the system is calculated using the time-evolution of the cluster population as input to the hydrodynamic description. (vi) The gelation time at which the gelation (fluid-solid) transition occurs is calculated analytically as the time at which the clusters connect into a random close packing and the low-shear viscosity diverges.

ii.2 Derivation

Master kinetic equation for aggregation with reversible bonds

Let us consider the master kinetic equation which governs the time evolution of the concentrations of clusters of any size present in the system as a result of the microscopic two-body association and dissociation processes:


where is the number concentration of aggregates with particles in each of them. is the rate of association between two aggregates, one with mass and the other with mass , while is the rate of dissociation of a aggregate into two aggregates and . The first term expresses the ”birth” of clusters with mass , the second expresses the ”death” of clusters with mass due to aggregation with another aggregate. The last two terms express the ”death” and ”birth” of -aggregates due to aggregate breakup, respectively. Instead of considering the two dissociation terms in the master equation explicitly, we can account for dissociation in an effective way by replacing the association constant with an effective size-independent rate constant and dropping the breakup terms in the master equation. If association is controlled by diffusion, as we are going to see in the next section, the rate of association is in good approximation independent of the sizes of the two colliding clusters. Further, we also assume that dissociation is also independent of the clusters. These simplifications are indeed justified if we assume that two clusters aggregate by forming a bond between two particles protruding on the respective surfaces, such that the association/dissociation kinetics between any two clusters can be effectively described by means of . The new master equation under these simplifications reads as:


Upon discrete-Laplace transforming this equation (26), the analytical solution for the time evolution of the cluster mass distribution (CMD) reads:


denotes the number per unit volume of monomer particles in the colloidal sol at . is the characteristic aggregation time or lag time and is equal to:


During this lag time, aggregation is slow because of bond breakage, and the formation of large clusters is unfavorable because with . Hence, for we can safely assume that the system is mainly composed of monomers and of doublets or colloidal dimers. If the bond dissociation process is stochastic, as we can anticipate, then after a time longer than the formation of stable bonds is possible. This has to be interpreted in a stochastic sense as the probability of large fluctuations around the average dissociation rate increases with time thus making possible the stochastic formation of long-lived bonds over a long time.

In the next section we derive an analytical expression for exploiting the fact that for only monomers and dimers are present in the system.

Effective association rate accounting for dissociation

We start by considering the kinetics of reversible association between two colloidal particles to form a dimer


The association rate be denoted by and the dissociation rate by . If we denote with the concentration of monomers at time and with the total concentration of monomers at , the evolution of is governed by:


where we made use of the conservation condition: , with the concentration of dimers. With the initial condition , Eq.(6) has the following solution:


with and . One should note that has dimensions of an inverse time, while has dimensions of because it is the rate constant of a bimolecular second-order reaction, whereas is the rate constant of a unimolecular or first-order reaction.

With only monomers and dimers as for , the temporal evolution of the average hydrodynamic radius of the system as measured by dynamic light scattering (DLS) is given by (27)


where is the intensity of light scattered by a population of pure monomers of radius , and is the intensity of light scattered by a population of pure dimers having a hydrodynamic radius . can be calculated by substituting Eq. (7) together with the conservation relation into Eq.(8). To first order in , the resulting expression reads as:


Upon taking the derivative and rearranging terms we obtain the standard form


The truncation to first-order in time implies that we are neglecting the equilibrium plateau that ultimately is reached according to the law of mass action. Rigorously, this approximation is valid for as discussed in the Appendix A. While keeping this in mind, it is instructive to consider its predictions also outside the rigorous regime of validity. By comparing the previous two expressions, we are now able to obtain the effective association rate accounting for bond-dissociation:


where . Since , at the particle length scale one has , while at the same time for spheres. Hence, in good approximation and to make our formulae more transparent, in our model we take , which is not going to change our results neither qualitatively nor quantitatively. With this result we can also identify the lag time:


This framework allows us to account for dissociation in the kinetics of aggregation between colloidal particles, in an effective way. It contains indeed both the association rate and the dissociation rate .

These rates can be estimated by solving the stationary equation of diffusion for the two particles in the frame of one of the two taken to be the origin. In the case of association, upon assuming stick-upon-contact as for short-range attraction one has: with the boundary conditions at and at . At steady-state, the solution for the rate of collision per unit volume or flux follows upon integration as: , which, upon using the Stokes-Einstein relation, leads to the Smoluchowski rate . The association rate is therefore independent of the size of the two particles or clusters that aggregate. This fact implies that larger particles (or clusters) aggregate at the same rate as smaller particles because the lowering of the diffusivity brought about by the larger size is exactly compensated by the increase in the collisional cross-section. When the attraction range cannot be neglected, one has to solve the diffusion equation for two particles in the field of force of a rectangular well of depth and width . The result is (28):


where is the mutual diffusion coefficient of the particles. For short-range potentials one recovers the Smoluchowski rate which we are going to use throughout this work. We should also mention that hydrodynamic interactions and elastic deformation effects of polymer-functionalized surfaces might play a role as well in the very short ranged limit. Since we cannot accurately model the latter effect we choose here to use the classic Smoluchowski rate theory where the slowing down of the rate brought about by hydrodynamics near the surface cancels, approximately, with the speeding up brought about by the finiteness of the attraction range. For a discussion of this effect see Ref.(29).

In a similar fashion, the dissociation rate can be estimated by solving the steady diffusion equation for the two bonded particles in the field of the attraction potential. At steady-state one calculates the (Kramers) rate of escape of one particle from the attraction well which drives the dissociation event. What changes with respect to the association process is obviously the boundary condition, since now it is assumed that the two particles start in a quasi-equilibrium steady state in the attraction well. In good approximation, the result for the average dissociation rate is:


In Ref.(28) a theoretical justification for this formula and its derivation can be found. Here we chose this simplified form and neglect non-essential prefactors in order to be consistent with our previous characterization of the colloidal systems under study (25). Eq.(14) is not accurate, however, as soon as one deals with shallow attraction energies . Even in the limit this formula predicts a finite dissociation rate although it is clear that no bond can be present to start with because the particles are hard-spheres and dissociation is therefore instantaneous as the spheres collide. This artifact is due to the approximation under which Eq.(14) is derived, whereby Kramers assumed in its derivation that the attraction well must be significant in order for the two particles to be in a quasi-equilibrium steady-state in the well (30). To overcome this problem rigorously, one should solve the full time-dependent diffusion equation which however would undermine the analyticity of our approach. Hence here we propose the following semi-empirical formula which interpolates between the Kramers formula for and the limit at :


In our calculations below, we are going to use and , With this choice, the dissociation rate is equal to the Kramers formula for all attractions down to and below it rapidly increases and diverges at . This interpolation formula is plotted in Fig.2 together with the Kramers formula for comparison.

Figure 2: (color online). Comparison between the Kramers dissociation rate (continuous line) given by Eq.14 and the interpolation formula Eq. 15 which interpolates between the Kramers formula and the limit.

With these identifications, the expression of the lag time, with the explicit dependence on , is given by:


The lag time is thus a function of the competition between microscopic association and dissociation kinetics. In the limit controlled by association , the lag time is set by the time of diffusive transport as for irreversible diffusion-limited aggregation: . In the opposite limit where dissociation is controlling, , the lag time goes as . Finally, when the condition , is exactly satisfied, which fixes , such that the initial individual frequencies of the forward process (association) and of the backward one (dissociation) are equal, the lag time scales as:


and it is inversely proportional to through the inverse of the Smoluchowski aggregation rate. This is a physically meaningful outcome because in this regime an increase of causes the speeding up of the diffusive transport which reduces the lag time. The physical meaning of this result is that, in the regime of equilibrium aggregation, the kinetics is controlled by the activated stochastic jump of the particles out of the attractive well which is the kinetically limiting process. Upon reducing the attraction, the lag time increases because the formation of bonds requires stochastically a longer time. Viceversa, upon increasing the attraction the lag time gets reduced because it becomes stochastically more likely to form long lived bonds on shorter times. Although it is tempting to identify the condition with the condition of detailed balance straightaway, this is not a rigorous identification. As discussed in the Appendix A where we refer to Tolman’s definition (31), this condition rather corresponds to the microscopic reversibility in the early time limit.

Finally, it should be noted that, once that , the aggregation kinetics enters an extremely fast regime which is independent of the lag time and hence is the same for all (finite) attractions. During this fast regime the kinetics is relatively insensitive to the microscopic details of cluster aggregation and it is very difficult to detect the effect of the microscopic details in the macroscopic properties which evolve very rapidly towards the solid state.

Analytical solution for the clustering kinetics with reversible bonds

We now have a connection between the mesoscopic time-evolution of the clusters and the microscopic interactions between colloidal particles which can also be allowed to vary with time. The CMD can be used to derive the time-evolution of quantities such as the average cluster size and the effective volume fraction in the system defined as the fraction of volume occupied by the clusters at time . In general, the radius of the sphere enclosing the cluster is given by: where is the cluster fractal dimension or mass-scaling exponent: , where is a prefactor of order unity. Using Eq. (3) for , the average cluster radius can be readily calculated:


Here, denotes the polylogarithm (32) of order of the variable . Closed form expressions for the polylogarithm of negative order are known only for the special case which corresponds to . Using this form we have that , which is a significant simplification. Extrapolating this result to higher values of we obtain a good closed-form approximation:


One can check that this approximation is reasonable quantitatively and it is always qualitatively correct in the whole range of and . This formula implies that in the aggregation process there is a lag time of order , during which aggregation events are stochastically rare because of the low rate of successful collisions leading to bond-formation, where successful collisions are those which do not result in immediate thermally-activated bond rupture. If is a small number, meaning that the lag time is short as for irreversible aggregation, then the kinetics transitions after short time to the growth law that has been reported experimentally for irreversible colloidal aggregation in the past (33). In the limit of equilibrium aggregation where microscopic reversibility is satisfied (see section IV.a and the Appendices for the connection between thermodynamic equilibrium and ), by combining the above expression with Eq.(16), this treatment gives:


which in the asymptotic limit correctly recovers the well known Lifshitz-Slyozov (34) scaling for the growth rate in the coalescence (coarsening) regime of phase separation following nucleation under equilibrium conditions. The link between nucleation and phase separation is discussed more in detail in section IV.a. Hence Eq.(19) is important because it covers all limits of colloidal aggregation kinetics, from irreversible aggregation to nucleation at equilibrium, and provides theoretical justification to many experimental observations in the past.

Similarly, the effective cluster volume fraction is given by:


This definition of the effective cluster volume fraction is the most used in the colloidal gelation literature (35); (36). An alternative would be the mass-weighted average volume fraction which is obtained by introducing a factor in the sum in the numerator and is more used in the context of polymer gelation where the molecular weight is the key parameter in the gelation process. The mass-averaging is however less used with colloidal gelation as here the key parameter which controls the gelation process is the cluster size (and, of course, its cube which is the volume) rather than the mass. Using , where is the volume fraction of the colloidal gas and the volume of a single particle, and Eq.(3) for the CMD, we obtain:


Eq.(22) gives the effective cluster volume fraction as a function of the time-dependent interaction (accounting also for dissociation) embedded in the characteristic aggregation time . The polylogarithm of order with can be very accurately approximated as . Then the volume fraction occupied by clusters after some manipulation becomes:


Linking the clustering kinetics with the macroscopic viscosity

Consistent with our main approximation of treating the clusters as renormalized spheres occupying an effective rescaled volume fraction , we now describe the effective viscosity of the system as a function of . This treatment applies to fractal clusters as well in the regime where the screening of the hydrodynamic interactions from the interior makes them behave like effective spheres or spheroids, which is a known fact (38). The viscosity of the system can be estimated by treating the clusters as effective hard spheres since their hydrodynamic behavior is very close to that of hard spheres even for fractal clusters and the short-range attraction has little effect on the hydrodynamic viscous dissipation. Under these assumptions, it is possible to obtain analytical expressions in closed for for the viscosity over the entire range up to the cluster close packing where the system arrests. This approach has the advantage of providing analytical scaling laws for the gelation time as a function of the controlling parameters.

The differential effective medium theory allows us to calculate the viscosity of a dense suspension starting from Einstein’s method for calculating the viscosity of a dilute suspension of spheres. Because of the assumption of dilute and non-interacting particles, one first obtains a linear dependence, i.e. the Einstein formula which accounts for the hydrodynamic dissipation of a single cluster (39). Upon introducing a small increment of particles in the system and accounting for their mutual correlations, it is possible to account for the many-body hydrodynamic interactions as well as for excluded volume effects in an effective way, leading to the following expression(40):


This equation is a key result of this work. Here, is the random close packing fraction of spheres at which the viscosity becomes infinite. The latter value is the object of many detailed studies aiming at its precise definition. In particular, this value can vary depending on the size polydispersity (41) and on the particle interactions (42), in the range . For our scope these differences are irrelevant and we have checked that they do not minimally affect the qualitative predictions of our model. Hence, consistent with our hard-sphere approximation in the viscosity calculation, we take . Also, one should note that upon approaching the regime , the system undergoes a glassy dynamical arrest where the clusters become caged by their neighboring clusters. Within this regime the viscosity still increases, as a power-law of according to Mode-Coupling theories (43) and with an exponential dependence according to Adam-Gibbs theories (44), before diverging at the random close packing. These scenarios could be implemented in our frameworks, in future work, to provide a more detailed description in the glassy regime. Here we are interested in the overall kinetics of the process and defer this type of detailed analysis of the glassy regime to future work.

Upon replacing this expression in the above expression for the viscosity, we obtain a closed-form expression for the viscosity as a function of the microscopic interaction, of the volume fraction, and of time:




Predictions of Eq.(25) for the viscosity as a function of time are shown in Fig.3 for different values of the attraction energy. It is seen how the viscosity diverges at a well defined gelation time , which is the time at which the effective volume fraction occupied by clusters reaches the random close packing fraction . For shallow attractions, it is evident that the viscosity remains constant and not too much above the value of the initial colloidal gas for a long time before it starts to increase. This observation has very important implications. In most experimental studies in the past the separation between gelling and non-gelling systems has been defined somewhat arbitrarily. Indeed it is customary in the experimental practice to establish that a system does not gel under certain conditions if it remains in a fluid state and does not aggregate significantly over a chosen period of time. Clearly, in this way the choice of the time span might be such that the observation time is shorter than the lag time, i.e. and states that would gel after a time might be improperly classified as non-gelling. The expressions reported here can therefore be of help to experimentalists in establishing a more rigorous criterion for drawing experimental state diagrams of colloids.

Figure 3: (color online). Low shear viscosity calculated as a function of time with varying attraction energy according to Eq. 25 and using Eq. 23 for the cluster packing fraction for the case . The existence of a lag time for all attractions, after which the aggregation is very fast for all attractions, is evident. and . From right to left: , , , , .

Iii Gelation time

We define the gelation time as the time at which , where is the zero frequency viscosity. Of course, there are many criteria that have been proposed in order to identify the gelation point. One that is widely used is the Winter-Chambon rheological criterion based on the elastic moduli (22), which requires a thorough characterization of the frequency response. Since here we are interested in the low-frequency behavior and will present experiments done in this limit, the is a stringent criterion to identify the liquid to solid transition. In the present context, since we are focusing on the low-frequency behavior and we can deal with the viscosity in closed form, we choose the criterion which implies that at the gelation point the system is connected and does not flow at least over a long time scale.

From Eq.(25) it is possible to determine the gelation time as a function of the other microscopic parameters, such as the attraction, the aggregation rate and the volume fraction. By setting the content of the main bracket in Eq.(25) equal to zero, in theory one could solve for the gelation time for an arbitrary . In practice, the equation cannot be solved analytically for an arbitrary . By studying the solution for a few special cases such as and , and dropping high-order terms in , we obtain the following approximate formula:


Hence the dependence of is the same as for given in Eq. (16). As one could expect, also from the consideration of Fig. 3, the gelation time is proportional to the lag time . Using Eq. (12), we can relate the gelation time to the microscopic association and dissociation processes:


Here we wrote to make the full dependence on explicit. Let us discuss the various limits of this formula as a function of the attraction energy first. When the formed bonds are frozen in for a long time and . In this limit we get a gelation which depends only on the diffusive transport and on :


In the opposite limit of shallow attraction, , we get:


Finally, when the early-time microscopic reversibility condition is satisfied, we have that the gelation time follows an Arrhenius law as a function of the attraction: . This law implies that the gelation time, in this regime, increases upon increasing because the bonds become more short-lived at higher which slows down the aggregation process. Therefore, it is clear that the temperature affects the gelation in a very different way depending on the strength of the binding energy . In particular, in the regime of strong binding close to diffusion-limited aggregation, the dependence of the gelation time is governed by Eq. (17) and the gelation time decreases upon increasing because the diffusive transport is enhanced at higher . In the opposite limit of lower , instead, the gelation time obeys Arrhenius behavior and increases upon increasing because of the slowing down of aggregation caused by the enhanced thermal breakup of the bonds. Hence, both these predicted behaviors appear physically meaningful in the two opposite regimes.

Figure 4: (color online). Gelation time calculated as a function of the colloid fraction for different values of the attraction energy . . From top to bottom: , , , .

The gelation time as a function of the colloid fraction is shown in Fig.4 for . Upon increasing the attraction, the power-law decay with the exponent predicted in the limit of weak attraction gradually decreases and melds into the limiting power-law at high attraction.

Figure 5: (color online). Gelation time calculated as a function of the attraction energy for different values of the colloid fraction . . From top to bottom: , , , .

The behavior of the gelation time as a function of the attraction is plotted in Fig.5. Three different regimes can be identified. At low attraction, , the gelation time decays very rapidly from its asymptotic value at and melds into a second regime where the decay with the attraction is exponential. Clearly the first regime is dominated by the scaling , whereas the second regime is controlled by the scaling , according to Eq.(30). This predicted exponential scaling regime is confirmed by earlier experimental observations(45). In both these regimes the dissociation rate dominates over the association rate, i.e. . Upon further increasing the attraction, however, the dissociation rate becomes increasingly smaller in comparison with the association rate as and the exponential behavior flattens out into a plateau where the gelation time is independent of . This latter regime recovers the diffusion-limited irreversible aggregation which is characterized by permanent bonds since dissociation is now infinitely slower compared to association.

Iv State diagram of attractive colloidal matter

It is possible to summarize the model predictions in a state diagram of attractive colloids. In order to be consistent with previous studies, we introduce the effective temperature customarily defined as (16):


where is the colloid diameter. Clearly, low values of correspond to high attraction energies, and there is also a close relation with our dissociation rate, with low values of corresponding to low values of the dissociation rate . The state diagram is plotted in Fig.6. Continuous lines are calculated for specified values of and refer to nonequilibrium boundaries. In practice, every continuous line separates systems that undergo gelation on a time scale , at higher (i.e. to the right of the curve), from systems that undergo gelation at , at lower (i.e. to the left of the curve). All curves are plotted for since most experimental observations of gelation in this regime report values of close to this value, at least for . We have checked that changing in the range does not alter the results qualitatively. From these curves it is evident that at least in the lower half of the state diagram the gelation process is a strongly nonequilibrium process and the transition from sol to gel depends crucially on the time-scale of the experiment. In particular, if the time of the experiment is short compared to the gelation time , gelation cannot be observed and the system appears liquid-like and at most composed of freely diffusing clusters. On the other hand, if the time of observation is long compared to the gelation time, a transition from a liquid-like material into a solid-like one will appear. Hence, in light of our results, the observation of so called ”equilibrium” clusters in the absence of gelation in purely attractive colloids (46) might have been due to the time scale of the experiment being short compared to the theoretical time scale of gelation for those conditions (and in fact the attractions reported in Ref. (46) lie well in the lower regions of our diagram). The situation might be different, however, in the case of charged colloids where the electrostatic repulsion plays a major role giving rise to further effects (47) that are not considered in our analysis. The role of the time coordinate on the gelation transition has been neglected in many previous studies of colloidal gelation, both experimental and computational, despite being a key control parameter in all nonequilibrium transitions.

Figure 6: (color online). State diagram of attractive colloids. Continuous lines represent gelation lines (see text for explanation). From left to right: , , , , , , , , , . MR (dashed) line: this line represents the set of points on the diagram for which the early-time microscopic reversibility condition holds exactly. B (dotted) line: binodal line. S (dotted) line: spinodal line.

In the state diagram we have also plotted the early-time microscopic reversibility condition (DB line in Fig.6) represented by the equality between the rates of the microscopic forward and backward process, in this case association and dissociation, respectively. The early-time microscopic reversibility condition is expressed by the equality: . Regions of the state diagram lying much below this line are certainly away from thermodynamic equilibrium and in fact it is below the line that we observe the most striking dependence on time. Above this line, the gelation lines at different tend to become more closely spaced together until they almost merge together in the top part of the diagram. Interestingly, the region where all gelation lines (corresponding to gelation times separated by up to 10 orders of magnitude) practically merge, coincides with the binodal (liquid-liquid) phase-separation line calculated here using a simple mean-field Bragg-Williams approach. Hence, from the point of view of kinetic theory, gelation in the upper region of the diagram above the binodal line, is an extremely unlikely occurrence even over extremely long times. From the point of view of equilibrium thermodynamics, no gelation can happen above the binodal because here the free energy of the system is controlled by the unfavorable entropy of mixing which favors the sol state (that is obviously more ”mixed up” than any aggregated state) and makes macroscopic aggregation thermodynamically forbidden.

iv.1 The metastable region: nucleation and liquid-liquid phase separation

In Fig. 6 we have also plotted the spinodal line and the metastable region (delimited by the binodal and by the spinodal) deserves careful analysis. Although we have plotted our nonequilibrium gelation lines also in this region, they should be taken as purely indicative because in the metastable region, gelation is replaced by nucleation leading to liquid-liquid phase separation. Indeed, the metastable region appears to be centered upon the early-time microscopic reversibility line which is a necessary condition for the system to be close to the microscopic equilibrium between association and dissociation and for detailed balance to be satisfied (see also Zeldovich (48) for the detailed balance principle within the context of phase separation). Under these conditions, Eq.(1) leads straightforward to homogeneous nucleation, as shown in the Appendix. Nucleation leads to the formation of compact clusters which can be seen with a simple calculation.

Being close to the microscopic reversibility line where detailed balance can be satisfied (see the Appendix), we can write down the free energy for the formation of a nucleus or cluster. As in nucleation theory, the nucleus is treated as a macroscopic object which allows us to formulate the free energy of a single cluster. The Gibbs free energy contains two contributions. One is the volume enthalpy arising from the bonds that are formed: where denotes the number of particles in the cluster and the bond energy, as usual, and is the mean number of nearest neighbors. The other term is the energy spent to create the interface between the cluster and the solvent. For a fractal the interface is intrinscally discrete and the effective surface can be estimated as the surface occupied by the particles in the outmost shell of the cluster: , where we used that the number of particles in the outermost shell of a fractal cluster is equal to . In the limit one recovers for the cluster surface.

Therefore, the surface energy is equal to: , where is the surface tension. The free energy of the cluster is given by: . The cluster growth at equilibrium occurs along the path of minimal free energy. For all values of the parameters involved in the free energy, minimization of with respect to in gives . We did not consider the -dependence of . However, since is a monotonically increasing function of (because the average density of a cluster increases upon increasing at a fixed size), it is clear that accounting for its dependence on leaves this result unchanged.

iv.2 The unstable region: nonequilibrium aggregation

Hence we have shown that close to thermodynamic equilibrium the clusters are compact objects with which is in agreement with many experimental and simulation results presented in the past. As a consequence in the metastable region in between the binodal and the spinodal line the system is more likely to undergo aggregation into compact aggregates with the growth law derived in section II.b.3, and there is no competition with gelation.

It is possible, however, that a solid-like state is formed following spinodal decomposition if the volume fraction in the dense phase reaches the critical volume fraction for the attractive glass transition predicted by Mode-Coupling theories, according to a mechanism that has been recently discussed in several studies (49); (8); (12). Although this mechanism is certainly a good candidate to explain a fluid-solid transition in the proximity of the spinodal line, its application is however limited to the region of the phase diagram close to the microscopic reversibility line. For deeper quenches well below this line, we have and the time scale associated with the rearrangement of the bonds on the way to equilibrium (i.e. along the minimization of ) is thus longer than the time scale on which new particles join the cluster to form new bonds. The incoming particles stick irreversibly onto the particles in the outer layers of the clusters which leads to according to the well-known mechanism of diffusion-limited fractal aggregation.

Summarizing, the emerging picture is that sufficiently close to the microscopic reversibility line in Fig.6, the clusters can minimize their free energy which leads to . With this input, our analytical solution to the master kinetic equation for aggregation gives the growth law , in agreement with the coarsening kinetics of spinodally decomposing systems (34). Upon departing from the equilibrium conditions towards higher attractions, the relaxation time over which the cluster minimizes its free energy becomes long compared to the time scale of bond-formation due to the microscopic imbalance between association and dissociation. As a consequence, because the energy minimization cannot be completed on the gelation time scale and aggregation proceeds with a faster kinetics given by according to Eq.(19). Hence, under these conditions deep inside the unstable region in the state diagram, gelation cannot be driven by spinodal decomposition because the kinetics associated with the spinodal-like coarsening is slower than the kinetics associated with nonequilibrium fractal aggregation. As we are going to see below, this is already the case with a relatively mild attraction energy such as .

V Comparison with experiments

v.1 Complex viscosity

The model for the steady shear viscosity of the system as a function of the attraction presented in section II.b.4 can be used within a generalized hydrodynamics approach which bridges the gap between the hydrodynamic (small ) and the kinetic (large ) regimes and thus provides predictions for the rheological response in the whole frequency spectrum. According to generalized hydrodynamics, the constitutive relation is written in the following Maxwell form (50):


where and are the and components, respectively, of the microscopic displacement field, and is the shear stress. Upon Laplace-transforming the above equation we obtain the complex viscosity as:


where represents the instantaneous () shear modulus and is the Maxwell relaxation time.

Below the transition where the system is fluid, we use experimentally measured values of reported in previous studies  (51) for all calculations. As , is the one of the solidified system and can be evaluated as the affine contribution to the shear modulus of a disordered lattice of harmonically-bonded particles. The latter, consistent with our picture, are the clusters present in the system at the time at which  (24). In fact, whilst the low- shear modulus is strongly affected by nonaffinity (52), the high- one reduces to the affine shear modulus given by the following mesoscopic theory (52): . All the parameters in this formula can be evaluated as follows. The harmonic spring constant is , where is given by Eq.(35), and nm is the hydrophobic attraction range. The lattice constant is given by the average cluster diameter. This can be calculated using the CMD Eq.(3) evaluated at the time at which , and the result is . is the number density of clusters at the random close packing and is given by . Finally, the coordination number at the random close packing has to be (53) . Hence, there are no nontrivial adjustable parameters.

The observable quantity which is measured in the experiments is the modulus of the complex viscosity, defined as which then gives:


From this expression we can extract the most interesting limit which gives when , i.e. at the fluid-solid critical point.

Figure 7: (color online). Model calculation of the complex viscosity of thermosensitive colloids using the theory presented here [Eq.(34)]. See text for the estimate of . , , . The attraction energy is calculated using Eq. (31) and . As in the experiment, the temperature is a ramp function of time according to . From right to left: , , , , .

Eq.(34), together with the microscopic model for the cluster evolution as a function of time and attraction, can be used to study the rheological response of the system upon varying the microscopic parameters. Fig.7 displays the results of model calculations. The theory predicts a sharp rheological fluid-solid transition as a function of time upon varying the attraction provided that the final attraction strength is larger than a threshold, i.e. for the attraction range used here , typical of hydrophobic attraction (25). For there is no such sharp transition and gelation is a very slow process such that the final solid state may not be attained at all on the time scale of observation. The curves in Fig.7 refer to and we have checked that changing between two and three has basically no effect on the macroscopic transition, apart from shifting the solid plateau to somewhat higher values upon increasing , because at dynamical arrest clearly is smaller for larger and .

Figure 8: (color online). Comparison between the complex viscosity from the theoretical predictions of the model (lines) and the complex viscosity measured experimentally (symbols). Volume fraction: () , () . Solid line: theoretical prediction using the model presented here with .

v.2 Experimental system and interaction

The colloidal particles used here are the same ones characterized in detail in Ref. (25). They consist of a 52 radius solid polystyrene core onto which a polymeric network of crosslinked PNIPAM with -dependent thickness (nm at and nm at ) is affixed. All electrostatic interactions are fully screened by the addition of 5.10 potassium chloride. Attraction is induced by the hydrophobic effect as the particle shell upon collapsing becomes hydrophobic upon increasing . The bonds due to hydrophobic attraction become increasingly long-lived upon increasing because the rate of bond dissociation decreases with , according to Eq.(15). As for two-level systems, is a sigmoidal (Fermi-type) function of and goes from zero up to a plateau value across according to the following equation derived in Ref. (25):


where is the entropy change across the transition (25). The latter parameter is an experimental input which is fixed by the surface chemistry of the particles used in the experiments and it controls the rate of variation of with . The behavior of Eq.(35) for the experimental system used here is shown in Fig.1.

The complex viscosity as a function of was measured at in a stress-controlled rotational rheometer MCR 301 (Anton Paar) where is varied using Peltier elements with accuracy of . The modulus of the complex viscosity was measured at a strain upon varying at a rate of . The solid fractions of the dispersions investigated are in the range which corresponds to initial volume fractions occupied by the particles before the onset of the gelation .

v.3 Comparison

In the theoretical calculation we take in agreement with previous studies (54). The comparison between the theory and the experiments is shown in Fig.8. It is seen that the complex viscosity is constant up to which is the value at which the attractive interaction sets in (25). Within this regime, the response, according to our viscoelastic model is completely dominated by the dissipative part, . The theory is able to capture the very sharp jump of which grows by several orders of magnitude within a fraction of degree Kelvin shortly after the onset of attraction. Thus, given the fact that no nontrivial adjustable parameter has been induced, the agreement of theory and experiment can be regarded as excellent.

Vi Conclusion

In conclusion, we proposed a theoretical, fully analytical model that bridges the microscopic physics of colloidal interactions and physical bond-formation and dissociation, with the macroscopic rheology and the colloidal gelation transition. According to the model predictions, for shallow attraction energy , growth occurs with a lag time . The size grows as , asymptotically recovering the Lifshitz-Slyozov law in the equilibrium limit . The gelation time decays as a power law of the volume fraction with an exponent which is both attraction-dependent and -dependent. In the limit of weak attraction, , whilst is found in the infinite attraction limit. The theory also yields a state diagram of attractive colloids. In the metastable region, an equilibrium-like nucleation scenario can be recovered as long as detailed balance holds, at least approximately. Although it is possible that gelation is driven by spinodal decomposition close to the spinodal line, identifying a gelation boundary which is independent of time for deeper quenches seems unfeasible and nonequilibrium aggregation with is a faster process than spinodal coarsening with . Our results are supported by experiments on a model suspension of thermosensitive colloids in which is increased from zero in the same sample by raising . The comparison shows that gelation already for is a sharp nonequilibrium fluid-solid transition with , and culminates with a solid random packing of clusters with infinite viscosity. The solidification is irreversible unless one switches off the attraction by reverting , which is a unique feature of the experimental system under study.

We have also presented experimental observations of colloidal gelation in a system of thermosensitive colloidal particles where the attraction is varied from zero up to a maximum in the same system by simply varying . The complex viscosity as a function of can be quantitatively modeled with the theory presented here and exhibits a sharp gelation transition with . The comparison between theory and experiments also indicates that already with a mild attraction energy of short-range attractive colloids undergo fast nonequilibrium gelation with and with no detectable hallmarks of spinodal decomposition. Our analysis suggests that the latter mechanism is more likely to play a role and affect gelation for weaker attractions close to the spinodal line where the rates of two-particle association and dissociation are close in value and aggregation occurs with .

Clearly, the framework presented here opens the unprecedented possibility of engineering functional materials that turn into solid on a desired time-scale and that can be switched in a fully controlled and reversible way between fluid and solid states.

Let us consider one-step association and dissociation processes by which one particle joins a cluster and dissociates from a cluster, respectively. According to the definition given by Tolman (31), the principle of microscopic reversibility is satisfied when for every cluster the absolute frequencies of particle attachment (to a cluster) is equal to the absolute frequency of particle detachment (from a cluster): , here represents the number of monomers in the system. Clearly, in the early stage of the aggregation process, if the latter is sufficiently slow, can be taken as a constant and equal to the initial concentration of monomers . Then, upon neglecting cluster size dependencies, we recover what we called the early-time microscopic reversibility condition . When this condition is strongly violated and the association process to form dimers is initially much faster than the dimer dissociation process, , then the formation of trimers etc. becomes a faster process. As a consequence of the fast trimer formation, the equilibrium concentration of dimers given by where is given by Eq. (7), cannot be reached. Then, the condition for cannot be satisfied. In turn, this fact implies that also the detailed balance condition which has to be satisfied , cannot hold since it is already violated at least for . In view of this argument, it is clear that the condition that we encountered multiple times and plays a crucial role in the development of our theory, is tightly connected with the detailed balance principle and it may be regarded as a necessary (although probably not sufficient) condition for detailed balance to be satisfied. Hence, this argument also provides a simple but robust quantitative criterion to assess whether the aggregation process is nonequilibrium or not, with the inequality marking the crossover from the equilibrium aggregation (nucleation) dynamics into the nonequilibrium one.

Finally, we should interpret the approximation underlying Eq.(9) within this framework. The truncation to first-order in time is equivalent to neglecting the kinetic equilibrium plateau for the dimers and to restricting the analysis to the initial growth stage only which occurs before the plateau. The time at which the plateau is reached can be easily estimated from Eq.(7). Since the reaches a plateau when its argument is , this gives for the time required to establish the law-of-mass action kinetic equilibrium for the dimers. If the arrival frequency of the third particles is , when then only the linearly growing part of from Eqs.(7)-(8) is necessary to arrive at . Then the inequality to be satisfied for our approximation to be fully justified is which can be rewritten as . Hence, any result presented here in the regime should be regarded as merely illustrative. This fact is however not worrisome because in that regime the crossover into nucleation takes place and our kinetic theory anyways has to be replaced by nucleation and coarsening at some point.

In this Appendix we shall see how the classical homogeneous nucleation is recovered when the detailed balance condition is satisfied and the aggregation occurs with by minimizing the free energy of clustering according to the mechanism discussed in section IV.a. Recall that the master kinetic equation for aggregation and breakup reads as:


Let us assume that detailed balance is satisfied such that aggregation occurs under equilibrium conditions by minimizing the free energy . Minimization gives , mainly to maximize the mean coordination number . Then only detachments of individual colloidal particles are relevant whereas body-fragmentation plays no role because it involves the breakup of many bonds, which is a much slower process. Furthermore, clusters form as a result of spontaneous fluctuations and the growth of each cluster is independent of the behavior of the others. Clearly, under these conditions, the above equation can be rewritten accounting only for one-step association and dissociation processes:


where we have incorporated the constant factor inside the association rate constants, as the monomer concentration is assumed to be constant for the slow process. To shorten the notation we put etc. and rewrite this equation as:


Since the attraction is relatively weak and thermal dissociation is important, the principle of detailed balance is applicable in this limit. Hence, we now introduce the equilibrium or steady-state concentration of aggregates of size which is a Boltzmann function of the minimum available work (free energy) needed to form an aggregate of size : . Upon applying the principle of detailed balance we have:


These relations allow us to eliminate from Eq. (38) the quantities and and Eq. (38) becomes:


Let us now transform the discrete distribution into a continuous one where is a continuous variable expressing the cluster size. Denoting by the spacing along the axis between the neighboring sizes and , we have , etc. Then we get:


Since is constant and , and vary little within the length , one can do an expansion in power series of up to the first non-vanishing terms, under the assumption that varies slowly and hence . The Taylor expansion is done on the terms in and and is centered on :


We now expand the second term on the r.h.s.: