Open Access Article

This Open Access Article is licensed under a

Creative Commons Attribution 3.0 Unported Licence

Henry A.
Lockwood
,
Molly H.
Agar
and
Suzanne M.
Fielding

Department of Physics, Durham University, Science Laboratories, South Road, Durham DH1 3LE, UK

Received
28th November 2023
, Accepted 14th February 2024

First published on 22nd February 2024

Motivated by recent experiments studying the creep and breakup of a protein gel under stress, we introduce a simple mesoscopic model for the irreversible failure of gels and fibrous materials, and demonstrate it to capture much of the phenomenology seen experimentally. This includes a primary creep regime in which the shear rate decreases as a power law over several decades of time, a secondary crossover regime in which the shear rate attains a minimum, and a tertiary regime in which the shear rate increases dramatically up to a finite time singularity, signifying irreversible material failure. The model also captures a linear Monkman–Grant scaling of the failure time with the earlier time at which the shear rate attained its minimum, and a Basquin-like power law scaling of the failure time with imposed stress, as seen experimentally. The model furthermore predicts a slow accumulation of low levels of material damage during primary creep, followed by the growth of fractures leading to sudden material failure, as seen experimentally.

In this work, we focus on the irreversible failure of gels and fibrous materials, motivated in particular by recent experiments demonstrating that the behaviour of a protein gel under sustained loading is strikingly similar to that of brittle solids.^{5,18} In these experiments, an initial regime of slow creep over several decades of time t later gave way to catastrophic material failure at a time τ_{f}. In the primary creep regime, the shear strain γ increased up to O(1) in a sublinear way, γ ∼ t^{1−α} with α < 1, with the shear rate decreasing as ∼ t^{−α}, reminiscent of Andrade creep in solids.^{19,20} In a secondary crossover regime for times 0.1 ≲ t/τ_{f} ≲ 0.9, the shear rate attained a minimum _{min} at time τ_{min}. The shear rate then increased dramatically as ∼ (τ_{f} − t)^{−1} in a final tertiary regime, diverging in a finite time singularity at τ_{f}, signifying catastrophic material failure. Across all three regimes, the shear rate fit the single master relation

(1) |

The same experiments demonstrated the failure time to decrease with imposed stress as τ_{f} ∼ Σ^{−β}_{0}, reminiscent of the Basquin law of fatigue in solids,^{21} and to scale linearly with τ_{min}via a Monkman–Grant relation,^{22}τ_{f} = τ_{min}/c with c < 1. Importantly, this suggests that sudden failure at time τ_{f} might be forecast once a shear rate minimum has earlier been observed at τ_{min}.

By combining rheology with optical and ultrasonic imaging,^{23} these experiments also tracked the strain field within the gel, in tandem with the global strain γ. In the primary creep regime, no macroscopic strain localisation was seen, although the authors noted that this does not rule out local rearrangements on scales below those imaged. Indeed, separate recovery tests suggested low levels of irreversible material damage during primary creep. Fractures then developed during the secondary regime, and grew to cause failure in the tertiary regime.

Motivated by these experiments, we introduce a simple mesoscopic model for the irreversible breakup of gels. Simulations in the step-stress protocol show it to capture the experimental phenomenology just described: (i) primary creep with decreasing shear rate ∼ t^{−α}, (ii) a secondary cross-over regime in which the shear rate attains a minimum at time τ_{min}, (iii) a tertiary regime in which the shear rate increases dramatically as (τ_{f} − t)^{−1}, leading to (iv) catastrophic failure in a finite time singularity at time τ_{f}, (v) a linear Monkman–Grant scaling of τ_{f} with τ_{min}, (vi) a Basquin-like scaling of τ_{f} with Σ_{0}, (vii) accumulating low levels of damage during primary creep, followed by (viii) the growth of macroscopic fractures leading to catastrophic material failure.

Although we shall mostly compare our results with the particular experiments of ref. 5 other experiments^{18} report the same phenomenology, with quantitative values of the exponents α, β however being non-universal across materials. Therefore, our focus will be on robustly capturing all the qualitative behaviour just described, while also characterising how the values of exponents such as α, β, etc. depend on model parameters.

To model the much slower process of gelation during a “waiting time” −t_{w} < t < 0 prior to the switch-on of stress at t = 0, we must of course model the (re)formation of gel strands. During gelation, therefore, we include an additional model ingredient: that after any element locally yields, it reforms with local strain l = 0 and a new energy threshold drawn from a distribution ρ(E) ∼ exp(−E/x_{g}) in which x_{g} is a model parameter. As an initial condition at t = −t_{w}, we assume a distribution of thresholds ρ(E). For x < x_{g}, the model captures ageing, with dynamics that become slower over time as the system gels. We take element reformation to be infinitely fast compared with this slow process of gelation, with τ_{ref} = 0. This is oversimplified: a fuller description might consider a finite reformation time τ_{ref} consistently during gelation, creep and yielding. However, our focus here is not to model gelation itself, but the subsequent creep and yielding of a gel, however formed, under stress.

During gelation, with immediate element reformation after breakage via an assumed τ_{ref} = 0, compared with the much slower timescale of gelation, the model is the same as the soft glassy rheology model.^{25} During creep, the model differs from the SGR model in having irreversible element breakage via an assumed τ_{ref} = ∞, compared with the much faster timescale of yielding.

So far, we have described the model as though the strain rate (t) were the same for all elements. In fact, to allow for strain localisation during material failure, we consider an array of S streamlines stacked across the flow gradient direction y of a shear cell of gap width L_{y}, with M elastoplastic elements on each streamline.^{26} Imposing force balance then allows us to calculate the shear rate profile (y, t), which becomes highly heterogeneous as the sample yields. We do this via two different algorithms. The first considers a solvent stress of small viscosity η = 0.05, additive to the elastoplastic stress,^{27} giving a total stress Σ = σ + η. The results shown are converged to the limit η → 0, suited to the relevant experimental regime, η ≪ kτ_{0}. The second has η = 0.^{26} The η = 0 algorithm was needed only to obtain the scaling of shear rate with time in the final tertiary regime, as the gel finally breaks completely. To seed the formation of strain bands, we initialise our creep simulations with a small heterogeneity: t_{w} → t_{w}[1 + εcos(2πy/L_{y})]. Indeed, in physical practice, most rheometers seed slight heterogeneity due to device curvature. For a fuller discussion of other possible sources of heterogeneity that might seed band formation, see ref. 28.

We choose units in which the attempt time τ_{0} = 1, local modulus k = 1 and gap width L_{y} = 1. We also set x_{g} = 1, thereby setting typical local yield strains of order unity. Unless otherwise stated we use S = 10 streamlines with M = 100000 elements per streamline, a perturbation amplitude ε = 0.1, and numerical time-step dt = 0.01. Parameters explored are then the imposed stress Σ_{0}, the waiting time t_{w}, and the noise temperature x.

Alongside the global strain signal γ(t), we also plot in Fig. 1(c) the degree to which the strain field becomes heterogeneous, i.e., shear banded, as a function of the flow-gradient coordinate y across the shear cell. Specifically, we characterise this degree of strain localisation, i.e., of shear banding, by defining the quantity Δγ(t). This is the standard deviation in strain across the streamlines stacked across y at any time t, normalised by the globally average strain γ(t). This remains small during primary creep, consistent with the low levels of irreversible damage evidenced via strain recovery experiments.^{5} The level of strain localisation then increases dramatically as the material yields, as seen experimentally.^{5}

As indicated by the dashed line, we define the failure time τ_{f} as the time at which the global strain attains a value equal to unity in Fig. 1(a). (In any simulations with η = 0, we instead define τ_{f} as the time at which the strain rate actually diverges. We have checked that these definitions give the same scaling of τ_{f} with Σ_{0} and τ_{min}, with only a small fractional quantitative difference.) This is plotted vs. the imposed stress Σ_{0} in Fig. 2, for several values of the waiting time t_{w} at a fixed noise temperature x. Collapse of the data when normalised as τ_{f}/t_{w} indicates a scaling τ_{f} ∼ t_{w}, consistent with the aging predicted by our model during gelation. The log–log axes of panel (a) reveal the Basquin-like scaling τ_{f} ∼ Σ^{−β}_{0} at low Σ_{0}, as seen experimentally with β = 5.5 in ref. 5 and β = 2.0 in ref. 18. This scaling predicts that any stress, however small, will always eventually cause failure. This has important implications for the shelf life of gel-based products (which are subject to gravity at least), and contrasts notably with the physics of soft glassy materials such as dense emulsions, colloids and microgels in which the timescale of yielding diverges at a non-zero yield stress Σ_{y}, with indefinite creep for Σ_{0} < Σ_{y}.^{29–32} The log–lin axes in (b) reveal τ_{f} ∼ exp(−mΣ_{0}) at larger Σ_{0}.

A notable achievement of the experiments in ref. 5 was to fit the shear rate across all three time regimes – primary, secondary and tertiary – to the single master scaling form of eqn (1). To explore whether such a fit holds here, we plot the strain rate (t) in three different ways in Fig. 3. Panel (a) shows the normalised strain rate /_{min}vs. normalised time t/τ_{f} on log–log axes, to emphasise the primary creep regime. The dotted line shows the power ∼ t^{−α}, with α = x. Panel (b) shows the same quantity, now on lin–lin axes to emphasise the secondary crossover regime. As seen experimentally, this secondary regime is seen for times 0.1 ≲ t/τ_{f} ≲ 0.9, with the shear rate attaining a minimum during it at t = τ_{min}. Finally, panel (c) plots the normalised shear rate /_{min}vs. normalised time to failure, (τ_{f} − t)/τ_{f}. In this figure, for consistency with the counterpart figure in the experimental ref. 5, we have inverted the time axis so as to progress towards failure from left to right. In this way, the limit τ_{f} − t → 0 at the right hand side of the figure corresponds to the limit in which t → τ_{f} and the sample finally catastrophically fails. The +1.0 power seen experimentally is indicated by the dotted line. The dashed line in all three panels shows the single experimental master scaling of eqn (1). This can be seen to perform well in the primary and tertiary regimes, but significantly less well in the secondary crossover regime.

Fig. 3 Strain rate normalised by its minimum value _{min}. (a) Plotted versus time t normalised by the failure time τ_{f}, with log–log axes to emphasise the primary creep regime, ∼ t^{−α}. Dotted line shows the power α = x = 0.5. (b) Plotted likewise, with lin–lin axes to emphasise the secondary crossover regime. (c) Plotted as a function of time remaining to failure τ_{f} − t, normalised by τ_{f}, to emphasise the tertiary regime of final failure. Note that the time axis is visually inverted so as to move towards failure from left to right. Dotted line shows the power +1.0. Dashed line in each panel shows the fit to the master scaling form of eqn (1) with μ = 0.019 and λ = 0.26. Imposed stress Σ_{0} = 0.1, 0.2…1.0 in curves black, red, green, blue, yellow, brown, grey, violet, cyan and magenta. Noise temperature x = 0.5, waiting time t_{w} = 10^{9}. M = 10000. η = 0.0. |

Another notable achievement of ref. 5 was to establish a linear Monkman–Grant relationship between the time τ_{min} at which the shear rate attains its minimum in the secondary regime, with the later time τ_{f} of catastrophic material failure in the tertiary regime. Such a relation suggests that the approaching time of material failure might be forecast in a given experiment, once the earlier shear rate minimum has been observed. Fig. 4 shows the same relation to hold across simulations performed for a wide range of x, t_{w} and Σ_{0}. The dashed line shows the linear relationship τ_{f} = τ_{min}/c with the same prefactor c = 0.56 reported experimentally.

Finally in Fig. 5 we plot the displacement profile U(y) at several times during creep and yielding. During early creep, the displacement profile remains essentially linear, corresponding to a near-uniform strain profile γ(y) = ∂_{y}U(y), consistent with the lack of strain localisation in the imaging experiments in this regime.^{5} In contrast, as the material later yields, strong strain localisation arises, as seen experimentally.

Despite these successes, our approach has several shortcomings. First, we assumed for simplicity a fast gel strand (re)formation time τ_{ref} = 0.0 compared to the slow timescale of initial gelation, but an infinite τ_{ref} compared to the faster timescale of yielding under stress. Future work might consider a finite τ_{ref} consistently during gelation and yielding. Arguably, however, our approach does not accurately model the process of gelation itself: the dynamics for t < 0 might better be seen as leading to a sensible initial condition for the creep simulations. As noted above, setting τ_{ref} = 0 throughout, including during straining, recovers the original SGR model. This has been used in earlier work^{37} to model creep and yielding of soft glassy systems, including gels, that flow after yielding. Second, the experiments of ref. 5 suggested part of the primary power law creep to arise from reversible viscoelastic squeezing of fluid through the gel matrix. In contrast, in our model all creep arises via local plasticity. Indeed, the origin of power law creep in amorphous materials is controversial, having been variously attributed to plastic rearrangements^{16} and linear viscoelasticity.^{4,5} Finally, we have modelled spatial heterogeneity only in the gradient direction y, whereas significant heterogeneity was also observed experimentally in the vorticity direction z. Future modelling efforts should increase the dimensionality simulated.

- W. C. Poon, L. Starrs, S. Meeker, A. Moussaid, R. M. Evans, P. Pusey and M. Robins, Faraday Discuss., 1999, 112, 143–154 RSC.
- P. J. Skrzeszewska, J. Sprakel, F. A. de Wolf, R. Fokkink, M. A. Cohen Stuart and J. van der Gucht, Macromolecules, 2010, 43, 3542–3548 CrossRef CAS.
- D. Bonn, H. Kellay, M. Prochnow, K. Ben-Djemiaa and J. Meunier, Science, 1998, 280, 265–267 CrossRef PubMed.
- F. Gobeaux, E. Belamie, G. Mosser, P. Davidson and S. Asnacios, Soft Matter, 2010, 6, 3769–3777 RSC.
- M. Leocmach, C. Perge, T. Divoux and S. Manneville, Phys. Rev. Lett., 2014, 113, 038303 CrossRef CAS PubMed.
- V. Gopalakrishnan and C. Zukoski, J. Rheol., 2007, 51, 623–644 CrossRef CAS.
- T. Gibaud, D. Frelat and S. Manneville, Soft Matter, 2010, 6, 3482–3488 RSC.
- V. Grenard, T. Divoux, N. Taberlet and S. Manneville, Soft Matter, 2014, 10, 1555–1571 RSC.
- J. Sprakel, S. B. Lindström, T. E. Kodger and D. A. Weitz, Phys. Rev. Lett., 2011, 106, 248303 CrossRef.
- J. H. Cho and I. Bischofberger, Soft Matter, 2022, 18, 7612–7620 RSC.
- S. Aime, L. Ramos and L. Cipelletti, Proc. Natl. Acad. Sci. U. S. A., 2018, 115, 3587–3592 CrossRef CAS.
- E. Moghimi, A. B. Schofield and G. Petekidis, J. Phys.: Condens. Matter, 2021, 33, 284002 CrossRef CAS.
- B. J. Landrum, W. B. Russel and R. N. Zia, J. Rheol., 2016, 60, 783–807 CrossRef CAS.
- S. N. Karobi, T. L. Sun, T. Kurokawa, F. Luo, T. Nakajima, T. Nonoyama and J. P. Gong, Macromolecules, 2016, 49, 5630–5636 CrossRef CAS.
- A. Pommella, L. Cipelletti and L. Ramos, Phys. Rev. Lett., 2020, 125, 268006 CrossRef CAS PubMed.
- H. Nechad, A. Helmstetter, R. El Guerjouma and D. Sornette, Phys. Rev. Lett., 2005, 94, 045501 CrossRef CAS PubMed.
- J. Rosti, J. Koivisto, L. Laurson and M. J. Alava, Phys. Rev. Lett., 2010, 105, 100601 CrossRef.
- J. Bauland, M. Leocmach, M.-H. Famelart and T. Croguennec, Soft Matter, 2023, 19, 3562–3569 RSC.
- E. N. D. C. Andrade, Proc. R. Soc. Lond. Ser. A-Contain. Pap. Math. Phys. Character, 1910, 84, 1–12 CAS.
- M.-C. Miguel, A. Vespignani, M. Zaiser and S. Zapperi, Phys. Rev. Lett., 2002, 89, 165501 CrossRef.
- O. H. Basquin, Proc., Am. Soc. Test. Mater., 1910, 625–630 Search PubMed.
- F. C. Monkman, Proc., Am. Soc. Test. Mater., 1956, 593–620 Search PubMed.
- T. Gallot, C. Perge, V. Grenard, M.-A. Fardin, N. Taberlet and S. Manneville, Rev. Sci. Instrum., 2013, 84, 045107 CrossRef.
- A. Nicolas, E. E. Ferrero, K. Martens and J.-L. Barrat, Rev. Mod. Phys., 2018, 90, 045006 CrossRef CAS.
- P. Sollich, F. Lequeux, P. Hébraud and M. E. Cates, Phys. Rev. Lett., 1997, 78, 2020 CrossRef CAS.
- S. Fielding, M. Cates and P. Sollich, Soft Matter, 2009, 5, 2378–2382 RSC.
- H. J. Barlow, J. O. Cochran and S. M. Fielding, Phys. Rev. Lett., 2020, 125, 168003 CrossRef CAS PubMed.
- R. L. Moorcroft and S. M. Fielding, J. Rheol., 2014, 58, 103–147 CrossRef CAS.
- T. Sentjabrskaja, P. Chaudhuri, M. Hermes, W. Poon, J. Horbach, S. Egelhaaf and M. Laurati, Sci. Rep., 2015, 5, 11884 CrossRef CAS PubMed.
- P. Ballesta and G. Petekidis, Phys. Rev. E, 2016, 93, 042613 CrossRef CAS.
- M. Siebenbürger, M. Ballauff and T. Voigtmann, Phys. Rev. Lett., 2012, 108, 255701 CrossRef.
- T. Divoux, C. Barentin and S. Manneville, Soft Matter, 2011, 7, 8409–8418 RSC.
- S. Pradhan, A. Hansen and B. K. Chakrabarti, Rev. Mod. Phys., 2010, 82, 499 CrossRef.
- E. Jagla, Phys. Rev. E, 2011, 83, 046119 CrossRef CAS PubMed.
- F. Kun, H. Carmona, J. S. Andrade Jr and H. J. Herrmann, Phys. Rev. Lett., 2008, 100, 094301 CrossRef CAS PubMed.
- Z. Halász, Z. Danku and F. Kun, Phys. Rev. E, 2012, 85, 016116 CrossRef PubMed.
- R. L. Moorcroft and S. M. Fielding, Phys. Rev. Lett., 2013, 110, 086001 CrossRef PubMed.

This journal is © The Royal Society of Chemistry 2024 |