Aging of Thermoreversible Gel of Associating Polymers

Associating polymers are heteropolymers which contain segments that attract each other via non-covalent interactions such as hydrogen bonding, hydrophobic and dipolar interactions, etc. Strauss (1989); Chassenieux et al. (2011); Zhang et al. (2018). In solution, these attractive segments can form intramolecular associations with similar segments of the same polymer or intermolecular associations between different polymers, depending on polymer concentration. Theoretical studies show that at low concentrations, intramolecular associations are favored leading to the formation of loops Cates and Witten (1986) and flower-like micelles Borisov and Halperin (1995), whereas at high concentrations intermolecular associations dominate and result in the formation of a macroscopic thermoreversible network in solution (in other regions of the phase diagram gelation may be replaced by phase separation) Semenov et al. (1995); Winnik and Yektaf (1997); Semenov and Rubinstein (1998); Rubinstein and Semenov (1998); Rubinstein and Dobrynin (1999); Dobrynin (2004). These physical gels exhibit viscoelastic behavior and can be classified as liquidlike or solidlike depending on the stability of the interchain associations that act as physical crosslinks. When the lifetime of these effective crosslinks is relatively short (on experimentally relevant time scales) the network behaves as polymer fluid with $G^{\prime}\ll G^{\prime\prime}$ where $G^{\prime}$ and $G^{\prime\prime}$ are the storage and loss moduli, respectively. Conversely, in the limit of large crosslink lifetimes (again, compared to experimental time scales), these gels exhibit solidlike behavior ($G^{\prime}\gg G^{\prime\prime}$) Rubinstein and Colby (2003). The structure of the gel depends on the number of attractive segments, on the strength of the interaction between them, and on their distribution along the polymer backbone.

Many physical gels show aging behavior Gomez-Solano et al. (2013); Secchi et al. (2013). For example, in an experimental study of methylcellulose gel (a thermoreversible network formed by aggregation of hydrophobic chain segments at elevated temperature), very slow reorganization of the network following its formation by temperature quench was observed by multispeckle dynamic light scattering, and it was shown that the timescale of density variations on length scales comparable to the wavelength of light, increases with the age of the gel Schupper et al. (2008). In this paper we use computer simulations of a simple model of associating polymers in implicit solvent Osmanovic and Rabin (2018) to study the aging of a thermoreversible gel. Similarly to ref. Semenov et al. (1995) where networks made of telechelic polymers were studied theoretically using mean field methods, we find that such gels are crosslinked via the formation of clusters of attractive segments (stickers), connected by soluble polymer segments. When the attractive interactions between stickers are sufficiently strong, we observe a dramatic increase of the lifetime of the crosslinks with aging time. We find that the observed aging is the result of structural reorganization of the clusters whose size increases with time and then saturates. Somewhat counterintuitively, we find that this reorganization takes place at fixed average number of polymers per cluster and is the consequence of the fact that the number of intramolecular associations and therefore of loops, increases with time.

We consider a system of $M$ polymers, of $N=50$ segments (beads) each. The backbone of the chain is connected by finitely extensible nonlinear elastic (FENE) bonds, with the potential

$$U^{FENE}=-0.5KR_{0}^{2}ln\left[1-\left(\frac{r}{R_{0}}\right)^{2}\right]+4\epsilon\left[\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\right]+\epsilon,$$

(1)

where we take $K=30.0$, $R_{0}=1.5$, $\epsilon=1.0$, $\sigma=1.0$ and truncate the $2^{nd}$ term on the rhs of the above equation at $2^{1/6}\sigma$ (thus taking only the repulsive part of the Lennard-Jones potential). Beads which do not have a FENE bond between them interact via the Lennard-Jones potential

$$U^{LJ}_{ij}(r)=4\epsilon_{ij}\left[\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\right]\hskip 28.45274ptr<r^{cut}_{ij}.$$

(2)

Here $U^{LJ}_{ij}(r)$ is truncated and shifted to zero at cutoff distance $r^{cut}_{ij}$. Each polymer contains two types of beads designated as stickers and non-stickers respectively, such that $\epsilon_{ij}=\epsilon_{s}$ if $i^{th}$ and $j^{th}$ beads are stickers and $\epsilon_{ij}=\epsilon_{ns}$ if at least one of $i^{th}$ and $j^{th}$ beads is a non-sticker ($r^{cut}_{s}$ and $r^{cut}_{ns}$ are the corresponding cutoff distances). We use LAMMPS Plimpton (1995) (Large-scale Atomic/Molecular Massively Parallel Simulator) to carry out Langevin dynamics simulations in the NVT ensemble. The simulation is performed in a box of size $51\times 51\times 51$ in units of $\sigma$, using periodic boundary conditions (the box size was chosen to be slightly bigger than the stretched length of the polymer). The motion of each bead is given by the Langevin equation, neglecting hydrodynamic interactions

$$m\ddot{\bf{r}}_{i}(t)=-\frac{\partial U}{\partial{\bf{r}}_{i}}-\zeta\dot{\bf{r}}_{i}(t)+\eta_{i}(t),$$

(3)

where $U$ (sum over all $U_{ij}$), $\zeta$ and $\eta_{i}$ are the total potential energy, bead friction coefficient and random thermal force due to implicit solvent, respectively. The rms amplitude of the random noise is proportional to $(\zeta k_{B}T/\Delta t)^{1/2}$, where $k_{B}$, $T$ and $\Delta t$ are Boltzmann’s constant, temperature and integration time step, respectively. All the time scales in the following are expressed in Lennard-Jones time units $\tau_{LJ}=(m\sigma^{2}/\epsilon)^{1/2}=1$ (the mass $m$, particle diameter $\sigma$, interaction parameter $\epsilon$, and $k_{B}$ are all set to 1).

We set the integration time-step to be $\Delta t=0.005$, friction coefficient $\zeta=0.02$, and fixed the temperature at $T=1$. Following ref. Osmanovic and Rabin (2018) we took $4$ stickers per polymer which are separated by $9$ non-sticker beads along the chain (the two non-sticker end segments are of lengths $9$ and $10$, respectively). To obtain an initially uniform polymer solution, we placed all the polymers in an array inside the simulation box. We assigned $\epsilon_{ns}=\epsilon_{s}=1.0$ and $r^{cut}_{ns}=r^{cut}_{s}=1.12\sigma$ (since the potential is shifted to zero at the cutoff, this yields purely repulsive interactions between all beads) and equilibrated the polymer solution for sufficiently long time to reach equilibrium. We then changed the interaction parameter between stickers to $\epsilon_{s}=6.0$ and the corresponding cutoff distance to $r^{cut}_{s}=2.5\sigma$ and continued to monitor the system through the processes of gel formation and aging.

Since our goal is to study the aging of gels, we began by choosing the polymer concentration well-above the sol-gel transition (we fixed the temperature at $T=1$ and, therefore, polymer concentration is our tunable parameter). The extent of gelation can be quantified by the parameter $M_{LC}/M$, where $M_{LC}$ is the number of polymers in the largest connected aggregate in the system Osmanovic and Rabin (2018). In order to define connectivity, we defined an association (bond) between two stickers if they are within a range of $1.5\sigma$ from each other; two polymers are considered to be crosslinked if they have at least one bond between them. A fully connected system corresponds to $M_{LC}/M=1$ where all the polymers belong to a single aggregate, whereas $M_{LC}/M=0$ corresponds to a perfect sol state where none of the polymers in the solution are bound.

Starting from a random initial condition (obtained by equilibrating the system under good solvent conditions for the polymers), we performed simulations for different volume fractions, $\phi$, and computed $M_{LC}/M$ as function of time for each $\phi$. We found that the size of the largest cluster rapidly increases with time and eventually saturates at a plateau value of $M_{LC}/M$ (a plot of the steady state value of $M_{LC}/M$ vs $\phi$ is shown in fig. S1A in SI). Since we are interested in the case where almost all the polymers belong to the network, we chose the polymer volume fraction $\phi=0.096$ for which approximately $98\%$ of the polymers belong to the connected network. This volume fraction exceeds the overlap volume fraction of polymers in good solvent (with $\epsilon_{s}=\epsilon_{ns}=1$) which we determined to be $\phi^{\ast}=0.059$ (not shown). Once the attractive interaction between the stickers is switched on, the network forms very fast and $M_{LC}/M$ reaches saturation at $t\approx 100\tau_{LJ}$ (see fig. S1B in SI).

In order to study aging we continued our simulations till $t=5\times 10^{5}\tau_{LJ}=500k$ ($k\equiv 10^{3}\tau_{LJ}$). Figure 1 shows a snapshot (2d projection of a 3d configuration) of all the stickers in the system at time $t=500k$. The non-sticker beads are not shown in the figure. We observe that the stickers form compact clusters of different sizes where the size of a cluster is defined as the number of stickers in it. To quantify the change of cluster size with aging time, we computed the distribution of cluster sizes at three different aging times, $t=5k,50k$ and $500k$ (see fig. 2A). Even though $M_{LC}/M$ saturates already at $t\approx 0.2k$, strong shift of the distribution towards larger cluster sizes was observed as aging time increased from $t=5k$ to $50k$, followed by a more moderate upward shift of the distribution with aging from $t=50k$ to $t=500k$. Since the total number of clusters decreases with time (see fig. 3) we conclude that the distribution ages due to dissociation of small clusters and the incorporation of the dissociated stickers into larger clusters. Note that in the latter stages of aging (between $50k$ to $500k$) the shift occurs almost exclusively at the lower end of the distribution and the largest cluster sizes are nearly unaffected by it. This concurs with our intuitive expectation that the size of a cluster is limited from above by packing constraints that arise from excluded volume repulsion between the non-sticker chains that are attached to it.

In order to obtain further insight into the physical mechanism that leads to the observed aging of the cluster size distribution, we note that each bond between stickers in a cluster is either intramolecular (between stickers of the same polymer) or intermolecular (between stickers of different polymers) in origin. Thus, a cluster can grow in size either by attachment of stickers of polymers which not yet bound to it, or stickers of polymers that are already bound to that cluster (the latter mechanism leads to the formation of intramolecular loops). Figures 2B and 2C show the distribution of number of different polymers and number of loops per cluster measured at three different aging times. Surprisingly, we did not observed any appreciable changes in the distribution of total number of polymers per cluster with aging time whereas the loop distribution shows exactly the same increasing trend with aging as the distribution of the number of stickers per cluster. Since a polymer with $p$ loops that emanate from a cluster contributes $p+1$ stickers to the cluster, these results suggest that the overall growth in the size of the clusters is the result of progressive loop formation as stickers of polymers that connect larger and smaller clusters, are released in the process of dissociation of the smaller clusters and are reattached to the larger ones, thus further increasing their size.

Next we examined the effect of aging on the lifetime of clusters in the network. We define the lifetime of a cluster as the time at which the number of the intermolecular stickers that were originally present in it at time $t_{0}$, falls below 2 (note that since the distribution of clusters remains fixed in steady state, this means that clusters are being continuously renewed). To measure the average lifetime of a cluster we identified all the clusters in the network containing two or more connected polymers at a particular aging time $t_{0}$ and monitored their composition (at time intervals of $2\tau_{LJ}$) for a period of $\Delta t=5k$ after $t_{0}$. The fraction of clusters that were not completely dissolved was computed as a function of time $t$ (see Figure 4). Beginning at short aging time, $t_{0}=5k$, $8\%$ of the initially present clusters have dissolved during $\Delta t=5k$ and continued to decrease after this time ($12\%$ have dissolved by $\Delta t=10k$ - not shown). At longer aging time, $t_{0}=50k$, only $2\%$ of the initially present clusters have dissolved during $\Delta t=5k$ (less than $4\%$ have dissolved by $\Delta t=10k$). At the longest aging time ($t_{0}=500k$), only $1\%$ of the clusters have dissolved during $\Delta t=5k$. Recall that we did not observe any significant changes in the distribution of total number of polymers in the cluster with aging (fig. 2B). On the other hand, we defined the stability of a cluster in terms of change of the identity of the polymers attached to that cluster and found that this stability strongly depends on the aging time (fig. 4). We conclude that the observed increase of the size of the clusters with aging (see fig. 2A) and of their stability are the result of the increase of the number of loops (fig. 2C). Thus, even though the number of polymers in a cluster does not change with time, the number of intermolecular associations between them increases with aging (as each polymer contributes an increasing number of stickers to the cluster) and the binding energy that holds the cluster together increases as the result. Since the dissociation of a sticker from a cluster is a thermally activated process, we expect the detachment time to increase exponentially with the binding energy of the sticker in the cluster Semenov et al. (1995); therefore, the stability of clusters is expected to dramatically increase with aging, as observed.

The configurations of a polymer in a gel can be classified according to the topology of intramolecular associations between its stickers. Figure 5 shows all the possible configurations of an associating polymer with four stickers. These configurations range from a fully extended state with no intramolecular associations between its stickers (configuration 1), to a completely collapsed state (configuration 11) where the four stickers form a single aggregate. A histogram that shows the probabilities of occurrence of these configurations in the gel at different aging times is shown in fig. 6. As aging time increases the probability of configurations (1-4) with no loops or with a single loop decreases and the probability of configurations (6-8,11) having two and three loops increases. The main change in the histogram takes place between $5k$ and $50k$ and only small aging effects occur later on, a trend similar to that observed in figs. 2A and 2C and in fig. 4. Configurations 6 and 11 that have the largest intramolecular aggregates (3 and 4 stickers, respectively) become dominant at longer aging times ($50k$ and $500k$). Chain entropy and excluded volume also play a role in determining the relative probability of different configurations. For example, configuration 7 is suppressed compared to configuration 6 even though they both have an aggregate of 3 stickers, presumably because it is more difficult to accomodate an isolated monomer on a loop in a neighboring intermolecular cluster than a monomer which resides near a free end of a polymer.

The size and geometry of clusters of stickers that crosslink the polymer network is determined by the interplay between binding energy that tends to maximize the number of nearest neighbors of a sticker (e.g., by maximizing the ratio of bulk to surface stickers in the cluster) and by packing constraints due to excluded volume interactions between non-sticker chains attached to the stickers. In order to get some insight about the geometry of the clusters, we computed the potential energy $\epsilon_{P}$ of a sticker in a cluster by adding the energy of its interactions with all other stickers within cutoff distance ($2.5\sigma$). We carried out this calculation for all stickers in all clusters and obtained a distribution of sticker energies in the system, which was then time-averaged over a time window of $5k$ after each aging time. Figure 7 shows the distribution of potential energy of the stickers at three aging times. Each peak corresponds to a different number of nearest neighbors of the stickers: the peak at $\epsilon_{P}=0$ corresponds to isolated (non-aggregated) stickers. A second peak at $\epsilon_{P}\approx-6$ corresponds to stickers that have only one nearest neighbor and its small amplitude reveals that isolated pairs of stickers are rare. We did not observe a peak that corresponds to two nearest neighbors at $\epsilon_{P}\approx-12$. A plausible explanation is that this happens due to the contribution to the potential energy from interactions with second nearest neighbors that have a range of possible radial separations (a plot of the radial distribution is shown in fig S2 in SI) from the sticker, therefore, a range of interaction energies. When these interactions are added to the nearest neighbor potential energy contribution and averaged, the second peak is smeared out. To confirm this assertion we computed the average potential energy contribution due to nearest neighbors only (with a cutoff of $1.5\sigma$) and observed a peak corresponding to two nearest neighbors (see fig. S3 in SI). Note that averaging over second nearest neighbor contributions does not effect the peak that corresponds to stickers which have only one nearest neighbor because such stickers mostly form pairs which do not have any second nearest neighbors. The contribution of second nearest neighbors is much smaller than that of nearest neighbors for stickers that have more than two nearest neighbors, and well-resolved peaks are observed for 3, 4, 5 and 6 nearest neighbors ($\epsilon_{P}\approx-18,-24,-30$ and $-36$, respectively). A downward shift of the distribution at high potential energies (small number of nearest neighbors) and an upward shift of the distribution at low potential energies (large number of nearest neighbors), is observed at longer aging times. Since increasing the number of nearest neighbors increases the cohesive energy of the cluster and therefore increases its stability, this observation provides additional insights about the connection between the structural reorganization of the network and the enhanced stability of its crosslinks (clusters of stickers), during the aging of the gel.

So far we have explored the aging of a gel made of polymers with strongly associating stickers ($\epsilon_{s}=6$ in units of $kT$). How sensitive is the observed aging behavior to changes of the interaction strength $\epsilon_{s}$ between stickers? To explore this question we decreased the interaction parameter from $\epsilon_{s}=6$ to $\epsilon_{s}=5$ and performed simulations at the same density ($\phi=0.096$) and temperature ($T=1$) as in the $\epsilon_{s}=6$ case. We first computed the fraction of polymers in the gel and found $M_{LC}/M\approx 91\%$, a value that is somewhat lower than $98\%$ found for the $\epsilon_{s}=6$ system (see fig. S4 in SI). We did not observe appreciable changes of the distributions of number of stickers, polymers and loops in clusters with aging time in the $\epsilon_{s}=5$ system (see figure 8). Even though all the distributions were quite broad, they were peaked at 2 stickers, 2 polymers and zero loops per cluster, respectively, from which we concluded that the most probable “crosslinks” in the network involve a single association between two polymers.

The fact that no significant aging was observed in Fig. 8 suggests that unlike the $\epsilon_{s}=6$ case, for $\epsilon_{s}=5$ the clusters reach equilibrium on time scale comparable with the time it takes a connected network to form from solution ($M_{LC}/M$ saturates in $0.3k$ - see fig. S4 in SI). Next, we considered the kinetics of dissolution of clusters in the $\epsilon_{s}=5$ system. No aging is observed in fig. 9 that shows the fraction of intact clusters as a function of time and almost $100\%$ of the clusters dissociated within $\Delta t=5k$ time interval, independent of aging time (recall that only $1\%$ of the clusters dissociated within this time interval in the $\epsilon_{s}=6$ system upon aging to $500k$).

Another comparison of relaxation kinetics in the two systems is shown in Fig. 10 where we plotted the average number of stickers in a cluster $\langle n(t)\rangle$ as a function of time. We observed that while in the $\epsilon_{s}=6$ system relaxation to steady state ($\langle n\rangle\approx 9$) takes $100k$ Lennard-Jones time units, in the $\epsilon_{s}=5$ system steady state is achieved ($\langle n\rangle\leq 5$) already after $1k$. We also computed the distribution of different polymer configurations in the system and found that for $\epsilon_{s}=5$ case, the most probable configuration for all aging times is configuration 1 in which there are no loops and therefore no intramolecular associations (see fig. S5 in SI). Since the enhancement of intramolecular associations with time was shown to be the root cause of aging effects in the $\epsilon_{s}=6$ system, the observation that intramolecular associations are suppressed in the $\epsilon_{s}=5$ system concurs with the absence of aging in this system.

In order to understand how suppression of intramolecular associations in the $\epsilon_{s}=5$ system affects the structure of crosslinks (clusters), we calculated the distribution of the potential energy of a sticker $\epsilon_{P}$ normalized by the value of the interaction parameter, $\epsilon_{s}$. In the absence of next nearest neighbor interactions one expects to see a distribution with peaks at integers $0,-1,-2,-3,...$ that correspond to the number of nearest neighbors of a sticker in a cluster (fig. 11). For $\epsilon_{s}=5$ the highest amplitude peaks occur at low values of $\epsilon_{P}/\epsilon_{s}$ ($0$ and $-1$), whereas for $\epsilon_{s}=6$ the highest peaks occur approximately at $-4,-5$ and $-6$. We also computed $\langle\epsilon_{P}\rangle$ as a function of time for both $\epsilon_{s}=5$ and $\epsilon_{s}=6$ cases (see fig. S6 in SI). The magnitude of difference between the steady state values of the average potential energies of $\epsilon_{s}=5$ and $\epsilon_{s}=6$ system is $\langle\Delta\epsilon_{P}\rangle\approx 13$. Since the time for thermally activated dissociation of a sticker from a cluster increases exponentially with $\langle\epsilon_{P}\rangle$ Semenov et al. (1995) (actually this energy is reduced by excluded volume repulsion between the non-sticker segments of the polymer than emanate from the cluster), this explains the dramatic increase of the stability of clusters in the $\epsilon_{s}=6$ compared to the $\epsilon_{s}=5$ system (see movies M1 and M2 in SI where a single cluster of six stickers is followed during $\Delta t=5k$ after aging to $t=50k$, for $\epsilon_{s}=6$ and $\epsilon_{s}=5$ cases, respectively).

In this work we studied the aging of thermoreversible physical gels of polymers that contain strongly associating segments (stickers) connected by long soluble chain segments. These polyamphiphilic molecules assemble into a network of connected clusters that act as effective crosslinks. When such a network is rapidly formed from solution (e.g., by a temperature quench), initially it contains many unpaired stickers and a broad non-equilibrium distribution of clusters of different sizes. We find that this distribution relaxes with time, as small clusters disintegrate and their stickers are absorbed by larger ones. The total number of clusters decreases and their size increases with time and this process is accompanied by increase in the number of intramolecular associations/loops (fig. 12A). In order to understand how the process shown in fig. 12A affects the connectivity of the network, in fig. 12B we present the histogram of the fraction of polymers that take part in one cluster, and those that bridge between $2,3$ or $4$ different clusters. We observe that the most probable configuration is of a polymer that bridges between $2$ clusters and that aging decreases the number of polymers that bridge between $3$ or $4$ clusters and increases the number of polymers that bridge between $1$ or $2$ clusters. We therefore conclude that while aging increases the size of the clusters and therefore increases their stability, it reduces the bridging between the clusters and therefore decreases the connectivity of the network.

While a detailed study of the viscoelastic behavior of such thermoreversible gels is beyond the scope of this work we would like to make some qualitative comments about possible effects of aging. Since shear moduli of polymer gels increase with the density of crosslinks Rubinstein and Colby (2003) we expect the high frequency plateau modulus that controls the response of our thermoreversible gels to small deformations, to decrease with aging. Since the characteristic stress relaxation time is related to the the process of removing a sticker from the cluster, we expect the frequency at which the response of the gel to shear deformation changes from viscous to elastic, to decrease with aging. The nonlinear elastic response is expected to be controlled by the stability of crosslinks (clusters) and therefore the resistance to large deformations (and fracture) should increase with aging.

Note that in some respects the observed aging in our system resembles glassy behavior in solids - a system out of equilibrium that can only reorganize by overcoming increasingly higher energy barriers Palmer et al. (1984). In our case the barrier corresponds to the energy cost of moving a sticker from one cluster to another and, since the energy barrier increases with cluster size Semenov et al. (1995), the process slows down with time as the average size of the clusters increases. The analogy is incomplete since, unlike in glasses, in our system the process eventually terminates and further reorganization stops when clusters reach their equilibrium size for the given sticker pair interaction parameter ($\epsilon_{s}=6$), temperature ($T=1$) and polymer volume fraction ($\phi=0.096$). For weaker sticker interaction ($\epsilon_{s}=5$) relaxation to equilibrium is much faster and takes place simultaneously with the formation of the network. For larger interaction parameters we expect aging time to increase dramatically; equilibration would require much longer simulation times and the system may become jammed on accessible time scales.

We expect the results of this work to be relevant to synthetic physical gels in which the strength of associations may be tuned by changing temperature Schupper et al. (2008), pH Lauber et al. (2017), etc. They may also be relevant to biomolecular aggregates of intrinsically disordered proteins such as FG nucleoporins that can form liquid-like microdroplets which remain stable for hours and eventually solidify Celetti et al. (2020). Aging has also been observed in ameloid-like gels of nucleoporins in which initial aggregation due to hydrophobic interactions between FG repeats was followed by formation of rigid $\beta$ sheet structures that stabilize the network Adera et al. (2010). Note that in our model the intermolecular associations between the polymers become stronger due to intramolecular association between stickers of the same chain because now each polymer contributes several stickers to the cluster, therefore stabilizing it. We propose that in the case of intrinsically disordered proteins enhancement of intramolecular associations and stabilization of clusters during the process of aging corresponds to formation of $\alpha$ helices or $\beta$ sheets and subsequent strengthening of intermolecular aggregates due to strong intermolecular interactions between these elements of secondary structure.

Acknowledgement: YR would like to acknowledge helpful discussions and correspondence with Eduard Lemke. This work was supported by grants from the Israel Science Foundation and from the Israeli Centers for Research Excellence program of the Planning and Budgeting Committee.