# A unified state diagram for the yielding transition of soft colloids

Stefano Aime∗ Laboratoire Charles Coulomb (L2C), Université Montpellier, CNRS,
Montpellier, France  Present address: ESPCI, Paris, France
<EMAIL_ADDRESS>Domenico Truzzolillo  Laboratoire Charles Coulomb (L2C),
Université Montpellier, CNRS, Montpellier, France  David J. Pine  New York
University  Laurence Ramos  Laboratoire Charles Coulomb (L2C), Université
Montpellier, CNRS, Montpellier, France  Luca Cipelletti∗ Laboratoire Charles
Coulomb (L2C), Université Montpellier, CNRS, Montpellier, France  Institut
Universitaire de France<EMAIL_ADDRESS>
###### Abstract

Concentrated colloidal suspensions and emulsions are amorphous soft solids,
widespread in technological and industrial applications and studied as model
systems in physics and material sciences. They are easily fluidized by
applying a mechanical stress, undergoing a yielding transition that still
lacks a unified description. Here, we investigate yielding in three classes of
repulsive soft solids, using analytical and numerical modelling and
experiments probing the microscopic dynamics and mechanical response under
oscillatory shear. We find that at the microscopic level yielding consists in
a transition between two distinct dynamical states, which we rationalize by
proposing a lattice model with dynamical coupling between neighboring sites,
leading to a unified state diagram for yielding. Leveraging the analogy with
Wan der Waals’s phase diagram for real gases, we show that distance from a
critical point plays a major role in the emergence of first-order-like vs
second-order-like features in yielding, thereby reconciling previously
contrasting observations on the nature of the transition.

The yielding transition of soft glassy systems is of great relevance both in
technological and industrial applications and at a fundamental level [1].
Despite profound differences in their microscopic structural features,
yielding occurs with very similar macroscopic features in systems as diverse
as colloidal and nanoparticle suspensions [2], emulsions [3, 4, 5], star
polymers [6] and microgels [7]. This suggests the presence of an underlying
general framework, which has been addressed in recent experimental,
theoretical and numerical works [8, 9], leading to contrasting results.
Measurements of the macroscopic viscoelastic properties suggest that yielding
develops progressively as the system is driven far from the linear regime [2,
3, 6, 7, 10]. Various models such as the soft glassy rheology [8], the mode
coupling theory [11, 12] and fluidity [13, 14, 15] or on-lattice [16] models
reproduce the evolution of viscoelastic parameters across yielding. Recent
experiments and simulations probing microscopic quantities indicate that
yielding is associated with an increase of particle mobility [17, 18, 4, 19,
20, 21, 22, 5, 23], suggesting that it may be described as a dynamic
transition between a quiescent, solid-like state and a dynamically active,
fluid-like state, bearing analogies with equilibrium phase transitions, an
approach similar to that used in the past to describe other flow-induced
transitions [24]. Note, however, that this description does not take into
account the ultra-slow relaxations that typically occur in soft solids even at
rest [25, 26]. These works suggested contrasting scenarios for the yielding
transition. Some systems exhibit features typical of a first-order transition,
such as a discontinuous jump of the particle mobility [20, 4, 21, 22, 5] or of
structural symmetries [27], the coexistence of dynamically distinct states
[20], and hysteresis [28]. By contrast, in other cases yielding is described
as a rather continuous transition [29, 19], with features such as sluggish
dynamics [17, 18, 4, 19], enhanced fluctuations [4, 30] and growing length
scales [19] typically associated with a second-order transition. Thus, the
nature of the yielding transition remains elusive: there is a dearth of
experiments and modelling addressing the mechanical response and the
microscopic dynamics of a class of soft materials sufficiently diverse to
allow for a general description of yielding.

Here, we establish a unified view of the yielding transition of repulsive soft
colloids by combining experiments probing both microscopic and macroscopic
quantities with theoretical modelling and numerical simulations. We
investigate samples of three kinds: concentrated suspensions of microgel
particles (M) and charged silica nanoparticles (N), and dense emulsions (E)
(see Methods for details). All samples exhibit qualitatively similar behavior
in oscillatory shear tests at frequency $\omega$ and at variable strain
amplitude $\gamma_{0}$, as exemplified by Fig. 1a for microgels. For small
enough $\gamma_{0}$, $G^{\prime}$ and $G^{\prime\prime}$, the storage and loss
moduli, are independent of $\gamma_{0}$, $G^{\prime}>>G^{\prime\prime}$ and
the stress amplitude $\sigma$ grows linearly with $\gamma_{0}$, indicative of
a predominantly elastic, linear response. As $\gamma_{0}$ is increased, a
gradual transition to the nonlinear regime is observed: $G^{\prime}$ and
$G^{\prime\prime}$ deviate from their low-$\gamma_{0}$ behavior, with
$G^{\prime\prime}$ going through a maximum and eventually exceeding
$G^{\prime}$. Deviations from a purely harmonic response become non-
negligible, as shown by the growth of the normalized third harmonic amplitude
$I_{3}/I_{1}$ of the stress response. At the largest strain amplitudes,
$\sigma$ grows sublinearly with $\gamma_{0}$ and both moduli decay as power
laws: $G^{\prime}\propto\gamma_{0}^{-2\nu}$ and
$G^{\prime\prime}\propto\gamma_{0}^{-\nu}$, with a sample-dependent exponent
$\nu$ in the range 0.6-0.75, see Supplementary Table SI1.

Figure 1: Viscoelasticity and spontaneous dynamics of a dense microgel
suspension. a): Oscillatory rheology for a microgel suspension at effective
volume fraction $\varphi=1.5$ (sample M40s). Left axis: first-harmonic storage
($G^{\prime}$, black squares) and loss ($G^{\prime\prime}$, white circles)
moduli vs strain amplitude $\gamma_{0}$, at fixed
$\omega=0.157~{}\mathrm{rad~{}s}^{-1}$. Right axis: first-harmonic stress
amplitude ($\sigma$, gold crosses) and normalized third harmonic component of
the stress signal ($I_{3}/I_{1}$, brown pluses). Rheological quantities evolve
smoothly from solid-like behavior (small $\gamma_{0}$, blue shades) to fluid-
like response (red shades). b): Correlation functions measured for the same
system at rest, for scattering vectors $0.2~{}\mu m^{-1}\leq q\leq 4.5~{}\mu
m^{-1}$, increasing from brown to green shades. The correlation functions are
fitted with compressed exponential functions, with compression exponent
$\beta_{0}=1.65>1$ (lines). c): Relaxation rate $1/\tau_{0}$ vs scattering
vector. Same color code as in b). The line is a fit with slope 1, indicating
ballistic dynamics.

The range of strain amplitudes over which rheological quantities signal the
transition from solid-like to fluid-like behavior is quite broad, making it
difficult to determine the nature of the yielding transition [4, 10]. To gain
a deeper insight on yielding, we couple rheometry to measurements of the
microscopic dynamics, using dynamic light scattering or differential dynamic
microscopy (see Methods). Both methods quantify the dynamics via the intensity
correlation function $g_{2}(\tau)-1$, which decays from 1 to 0 as microscopic
displacements grow beyond a length scale $\pi/q\approx(0.1-1)~{}\mu\mathrm{m}$
set by the scattering vector $q$.

The spontaneous dynamics measured at rest are similar for all samples, and are
well described by a slow, compressed exponential relaxation:
$g_{2}^{(s)}(\tau)-1=\exp\left[-2\left(\tau/\tau_{0}\right)^{\beta_{0}}\right]$,
Fig.1b, with sample-dependent values of the spontaneous relaxation rate
$1/\tau_{0}$ and of the exponent $\beta_{0}>1$, see Supplementary Table SI1.
These dynamics are ballistic, as indicated, for the microgels, by the linear
dependence of $1/\tau_{0}$ with $q$, Fig.1c. Similar spontaneous dynamics have
been reported for many other jammed or glassy soft samples at rest, and are
attributed to the slow relaxation of quenched internal stresses [25]. To
investigate the microscopic dynamics under deformation, we apply an
oscillatory shear with angular frequency $\omega$ and measure $g_{2}-1$
stroboscopically, for $\tau$ values a multiple of the oscillation period
$2\pi/\omega$. The dynamics probed by this echo protocol [31] are only
sensitive to irreversible rearrangements, either spontaneous or induced by
shear. Figures 2a-c reveal striking similarities of the overall behavior of
the correlation functions across all samples. Under small strain amplitudes,
the dynamics are independent of $\gamma_{0}$, while they increasingly
accelerate with growing strain at larger $\gamma_{0}$. Concomitantly, the
shape of $g_{2}-1$ evolves from a steep compressed exponential decay at low
$\gamma_{0}$ to a stretched shape at large $\gamma_{0}$. Data at all strain
amplitudes are very well fitted by the following expression:

$g_{2}(\tau)-1=\left\\{\chi\exp\left[-\left(\Gamma_{s}\omega\tau\right)^{\beta_{s}}\right]+\left(1-\chi\right)\exp\left[-\left(\Gamma_{f}\omega\tau\right)^{\beta_{f}}\right]\right\\}^{2}\,,$
(1)

where $\beta_{f}$ and $\beta_{s}$ are (sample-dependent) constants, whereas
$\chi$, $\Gamma_{s}$ and $\Gamma_{f}$ vary with $\gamma_{0}$. The
dimensionless relaxation rates for the slow and fast relaxation modes,
normalized by the oscillation frequency $\omega$, are designated by
$\Gamma_{s}$ and $\Gamma_{f}$, respectively.

Figures 2d-f show the strain dependence of the normalized relaxation rates,
$\Gamma_{s,f}$, and of the slow mode amplitude $\chi$, obtained by fitting Eq.
1 to the correlation functions of Figs. 2a-c. Three regimes can be
distinguished. For small enough strain amplitudes, $g_{2}-1$ relaxes through a
single, slow compressed exponential mode ($\chi=1$), with a stretching
exponent $\beta_{s}\geq 1$ (see Supplementary Table SI1) and a strain-
independent relaxation rate $\Gamma_{s}$ close to that at rest,
$1/\omega\tau_{0}$. For the microgels, oscillatory tests at
$\omega=0.157~{}\mathrm{rad~{}s}^{-1}$ and
$\omega=3.14~{}\mathrm{rad~{}s}^{-1}$ indicate no dependence of the slow mode
with $\omega$ (in physical units), further confirming that the dynamics
observed in this regime are unaffected by shear and simply correspond to the
sample spontaneous relaxation. As $\gamma_{0}$ exceeds a threshold value,
correlation functions become strain-dependent. A second, faster mode,
characterized by a stretched exponential relaxation, adds to the spontaneous
relaxation mode, whose relative amplitude $\chi$ rapidly decays from 1 to 0
with increasing $\gamma_{0}$. Finally, for sufficiently large $\gamma_{0}$,
$\chi\approx 0$: the correlation functions are well fitted by a single
stretched exponential relaxation, with $\Gamma_{f}$ increasing as
$\gamma_{0}^{n}$, with a sample-dependent exponent $1<n<8$ (red symbols in
Figs 2d-f). In this regime, we check for microgels that the fast relaxation
rate, in physical units, is proportional to $\omega$, as expected in the case
of dynamics fully dominated by rearrangements induced by strain oscillations.
Moreover, we find that for the microgels and emulsions $\Gamma_{f}$ scales as
$q^{2}$ (see Supplementary Figs. SI8-SI9), the hallmark of diffusive motion,
as also reported recently in simulations [21] and in experiments on other
kinds of microgels under large shear strain [23]. This shear-induced diffusive
behavior at large $\gamma_{0}$ is analogous to the dynamics of equilibrated
dense colloidal suspensions at rest [32, 33], and contrasts with the ballistic
behavior at small strain or at rest. Our experiments thus show that at the
microscopic level yielding corresponds to a transition between ultraslow,
ballistic relaxations at small $\gamma_{0}$ (unaccounted for in previous
works) and fast, diffusive relaxations beyond yielding. In analogy to the
recently reported abrupt change of microscopic quantities such as the particle
mean squared displacement or diffusivity [4, 21, 23], the amount of
irreversible rearrangements [5, 17, 20], and the size of avalanches [22], the
correlation functions measured in our experiments exhibit a marked change in a
narrow interval of $\gamma_{0}$, indicative of a transition sharper than for
rheological quantities, compare the stars and the vertical lines in Figs.
2d-f.

Figure 2: Yielding of soft solids as a dynamic transition. a)-c): Intensity
correlation functions under oscillatory shear, for concentrated microgels (a,
M2s), nanoparticles (b, N45%), and emulsions (c, E74%), plotted vs the
normalized time delay $\omega\tau$, with $2\pi/\omega$ the period of the
oscillations. $\omega=\pi~{}\mathrm{rad/s}$ (resp., $2\pi~{}\mathrm{rad/s}$)
for samples M2s and N45% (resp., E74%). The strain amplitude $\gamma_{0}$
increases from blue to red shades, spanning the rheological yielding
transition (see Methods for sample details and Supplementary Table SI1 for the
probed deformation ranges). Symbols: experimental data. Lines: fits using Eq.
1. The compressed exponent $\beta_{s}$ of the slow mode is 1.3, 1.9, 1.4 for
microgels, nanoparticles, and emulsions, respectively. The stretching exponent
$\beta_{f}$ of the fast mode is a fit parameter shared between all data for a
given sample, yielding $\beta_{f}=0.4$, 0.9, and 0.3 for M2s, N45% and E74%,
respectively. The fitting parameters for the same samples are shown in d)-f).
Left axes: $\gamma_{0}$ dependence of the normalized rates $\Gamma_{s}$ (slow
mode, blue squares) and $\Gamma_{f}$ (fast mode, red circles). Right axes:
relative amplitude $\chi$ of the slow mode (stars). Lines: numerical results
of the general model for the yielding transition discussed in the text, see
Supplementary Table SI2 for details on the parameter values. In d)-f), the
vertical lines indicate, from left to right, the onset of the increase of
$G"$, the onset of deviations from linearity of $\sigma(\gamma_{0})$, the
maximum of $G"$, the crossover between $G^{\prime}$ and $G"$.

To rationalize these findings, we introduce a simple model. The sample is
coarse-grained on a lattice; each lattice site is attributed a relaxation rate
$\Gamma_{i}$ that depends on both the spontaneous relaxation at rest
$1/\omega\tau_{0}$ and shear-induced rearrangements $\Gamma_{sh,i}$:

$\displaystyle{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}}$
$\displaystyle{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}=1/\omega\tau_{0}+\Gamma_{sh,i}}}$
(2a)
$\displaystyle{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{sh,i}}}$
$\displaystyle{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}=\frac{K}{\gamma_{0}^{-n}+N^{-1}\sum_{j}\frac{\alpha_{ij}}{\Gamma_{i}\Gamma_{j}}}\,,}}$
(2b)

where $K$ is a constant whose physical meaning will be discussed later, and
where the sum in the r.h.s. of Eq. 2b runs over the $N$ nearest neighbors of
site $i$, with $\alpha_{i,j}$ coupling constants between the dynamics of sites
$i$ and $j$. Equation 2a states that the overall relaxation rate is the sum of
two independent contributions: $1/\omega\tau_{0}$, the spontaneous relaxation
rate, and $\Gamma_{sh,i}(\gamma_{0})$, the site- and strain amplitude-
dependent rate of the additional relaxation induced by shear. A similar
additive rule has been invoked in mode coupling-based models [34, 35, 36],
which postulated $\Gamma_{sh,i}\propto\gamma_{0}^{n}$. However, this form of
the shear-induced relaxation rate yields a smooth growth of $\Gamma_{i}$ with
$\gamma_{0}$, rather than a well-defined transition. Instead, we propose in
Eq. 2b an alternative ansatz for the shear-induced relaxation rate. It is the
simplest expression that accounts for the following physical ingredients: i)
$\Gamma_{sh,i}$ should vanish for small strain amplitudes, implying
$\Gamma_{i}\approx 1/\omega\tau_{0}$ in the $\gamma_{0}\rightarrow 0$ limit;
ii) in the opposite limit of large $\gamma_{0}$, the dynamics should be
dominated by the externally imposed shear, implying $\Gamma_{i}\rightarrow
K\gamma_{0}^{n}$, as measured in our experiments for the fast mode; iii) in
the intermediate regime, the shear-induced dynamics should be ruled not only
by the external drive, but also by the interactions between neighboring sites,
which we expect to slow down the system relaxation, as modelled by the sum
term in the r.h.s. of Eq. 2b. The latter is chosen in the spirit of dynamic
facilitation models for the spontaneous relaxation of glassy systems [37],
where sites with a higher-than-average relaxation rate facilitate the
relaxation of neighboring sites.

We start by considering the mean-field version of the model, where the
coupling constants and thus the relaxation rates are taken to be identical for
all sites, $\alpha_{ij}\equiv\alpha$ and $\Gamma_{i}\equiv\Gamma$. The mean
field model can be solved analytically by recasting Eqs. 2 as

${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\left(\Gamma-\frac{1}{\omega\tau_{0}}\right)\left(\gamma_{0}^{-n}+\frac{\alpha}{\Gamma^{2}}\right)=K\,,}}$
(3)

with $\alpha$ an average coupling constant. This equation is formally
identical to the Van der Waals (VdW) equation of state ruling the vapor-liquid
transition of real gases, with pressure $p$ volume $V$ and temperature $T$ in
the VdW’s equation replaced by $\gamma_{0}^{-n}$, $\Gamma$, and $K$,
respectively. The spontaneous non-dimensional rate $1/\omega\tau_{0}$ and the
coupling constant $\alpha$ play the role of the molecular volume and molecular
interaction parameter in VdW’s law, respectively.

In experiments, the strain amplitude is the control parameter, typically
plotted on the $x$ axis. In Fig. 3a, we rather choose $\Gamma$ as the
abscissa, to emphasize the analogy of ‘iso-$K$’ solutions of Eq. 3 with VdW
isotherms in a $pV$ diagram. We find that in the mean field model $K$ plays a
key role in differentiating samples that exhibit a yield transition from
samples that are predominantly fluid-like at all $\gamma_{0}$, as illustrated
by the three curves of Fig. 3a. For $K$ larger than a critical value $K_{c}$,
solid line in Fig. 3a, $\gamma_{0}^{-n}$ decreases monotonically with
increasing $\Gamma$. This corresponds to the smooth growth — with no yielding
transition— of the relaxation rate of concentrated yet equilibrated colloidal
fluids upon applying a mechanical drive [38] . For $K<K_{c}$, by contrast,
$\gamma_{0}^{-n}(\Gamma)$ becomes non-monotonic (dashed line in Fig. 3a):
within a finite range of strain amplitudes, a unique value of $\gamma_{0}$ is
now associated with multiple values of $\Gamma$. In a VdW fluid, this feature
is associated with the vapor-liquid phase transition: upon compression, the
system jumps from the vapor branch to the fluid branch of the isotherm line at
a pressure set by the minimization of the free energy and corresponding to
Maxwell’s equal area rule. In our model, non-monotonic iso-$K$ curves are
associated with yielding. Starting from a solid at rest and increasing
progressively $\gamma_{0}$, the system descends the solid-like (left) branch
of the equation of state, corresponding to small and nearly constant $\Gamma$.
In the representation of Fig. 3a, portions of the equation of state with
positive slope are nonphysical, because they correspond to faster relaxation
rates attained at lower strain amplitudes. Thus, the system has to jump from
the left branch to the right, fluid branch, which constitutes yielding in our
model. We find that in the mean field limit of the model, the jump occurs at
the minimum of the iso-$K$ line, from point 1 to point 2 in Fig. 3a.
Introducing disorder smears the transition and, in the limit of small
disorder, the yield strain is shifted to smaller values, approaching a value
set by the equivalent of Maxwell’s equal area rule [39] (points 1’  and 2’ ).
Finally, the dotted line of Fig. 3a represents the critical iso-$K$: in
analogy to the VdW’s critical isotherm, it has an inflection point but no
local minimum. Here, it separates systems that are fluid-like at all
$\gamma_{0}$ from systems that are solid-like at small enough $\gamma_{0}$.

Figure 3: Theoretical and numerical modelling of the yielding transition. a)
Lines: mean field model, Eq. 3, for three $K$ values: $K=4~{}10^{-3}$ (sub-
critical, dashed line), $K=K_{c}=4.2~{}10^{-3}$ (critical iso-$K$, dotted
line) and $K=5.9~{}10^{-3}$ (super-critical, continuous line). The circled
numbers, arrows and horizontal lines indicate the yielding transition, see
text for details. Color shades and color bar: probability distribution of the
local relaxation rates $\Gamma_{i}$ in simulations of the model with disorder
and parameters chosen to reproduce the data for the microgels M2s shown in
Fig. 2d: $K=4~{}10^{-3}$, $1/\tau_{0}=2~{}10^{-5}~{}\mathrm{s}^{-1}$,
$\bar{\alpha}=9~{}10^{-8}$, variance of $\alpha_{i,j}=8.9~{}10^{-16}$, $n=3$,
$\omega=\pi~{}\mathrm{rad~{}s}^{-1}$. b): Probability distribution of the
local relaxation rate for strain amplitudes corresponding to the fluid (blue),
coexistence (gray) and solid-like (red) regimes, same colors as the
corresponding horizontal lines in a). c): Snapshot of the simulated system
corresponding to the bimodal distribution in b), showing the coexistence of
regions with low (dark shades) and high (light shades) relaxation rates. In b)
and c), the parameters are the same as for the model with disorder in a).

The mean field model, Eq. 3, describes yielding as a first-order transition
between two dynamically distinct states, accounting for both the linear and
the fully fluidized regimes. However, it fails to properly capture the gradual
onset of the fast-relaxation mode and the regime of intermediate strain
amplitudes where both modes coexist. Quenched disorder is known to smear out
first-order transitions [40, 41]. To explore the role of disorder in our case,
we solve numerically the full model, Eq. 2, using model parameters that fit
the microgels data of Fig. 2d (details in Methods). In the presence of
disorder, $\Gamma_{i}$ varies from site to site: representative probability
distributions
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
are reported for three strain amplitudes in Fig. 3b. In agreement with
experiments, three different regimes are seen: (1) under small strain
amplitudes,
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
is unimodal, peaked around a small value comparable to the relaxation rate at
rest; (2) under intermediate strain amplitudes,
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
becomes bimodal as a consequence of the appearance of a second, faster mode
characterized by a rate
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{f}}}$,
typically well separated from
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{s}}}$
and growing with $\gamma_{0}$; (3) under large strain amplitudes, the slow
mode vanishes and
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
is again unimodal, but is now sharply peaked around
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{f}}}\propto\gamma_{0}^{n}$.

We associate the bimodal nature of
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
at intermediate $\gamma_{0}$ with the coexistence of slow and fast relaxation
modes observed experimentally, which smears the transition with respect to the
mean field prediction (compare the dashed line and the distribution of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}}$
indicated by the color shades in Fig. 3a). A spatial map of the local
relaxation rates reveals that slow and fast relaxing sites form a coarse
structure (Fig. 3c), consistent with the spatial localization of highly mobile
regions observed in the single-cycle dynamics of sheared emulsions [4]. The
separation between the two modes allows one to extract from
$\rho({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{i}}})$
two well-defined values of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{s}}}$
and
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{f}}}$,
as well as the relative weight $\chi$ of the slow mode. A suitable choice of
the model parameters $K$, $n$, and of a log-normal distribution of $\alpha$
(see Supplementary Table SI2) reproduces the experimental strain dependence of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{s}}}$,
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{f}}}$
and $\chi$ (lines in Fig. 2d-f). The good agreement between experimental data
and numerical results over up to two decades in applied strain supports the
model and highlights that disorder is indeed at the origin of the dynamic
coexistence spanning a finite range of $\gamma_{0}$.

One of the most powerful consequences of VdW’s theory is the law of
corresponding states, predicting identical properties for distinct fluids,
provided that they all have the same pressure, volume, and temperature
relative to the corresponding values at the critical point. Inspired by the
law of corresponding states, we re-express Eq. 3 using reduced variables:

$\left({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}}}-\frac{1}{3}\right)\left({\gamma_{r}}^{-n}+\frac{3}{{{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}}}}^{2}}\right)=\frac{8}{3}K_{r}\,,$
(4)

where
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}=\Gamma/\Gamma_{c}}}$,
$\gamma_{r}=\gamma_{0}/\gamma_{0,c}$, $K_{r}=K/K_{c}$, and where the values of
the various parameters at the critical point, designated by the subscript $c$,
are given in Table 1. For the mean-field model, the coordinates of the
critical point are derived by imposing that both the first and the second
derivative of
$\gamma_{0}({{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma}})$
vanish, in analogy to VdW’s equation of state. For the model with disorder, we
use reduced variables obtained from Table 1 with the mean-field $\alpha$
replaced by the average value $\bar{\alpha}$ of the site-dependent
$\alpha_{i,j}$.

$\gamma_{0,c}$ | ${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{c}}}$ | $K_{c}$
---|---|---
$\left(\frac{27}{\alpha{{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\omega^{2}\tau_{0}^{2}}}}\right)^{1/n}$ | ${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\frac{3}{\omega\tau_{0}}}}$ | ${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\frac{8\alpha\omega\tau_{0}}{27}}}$
Table 1: Values of the strain amplitude $\gamma_{0}$, relaxation rate
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma}}$
and $K$ parameter at the critical point as predicted by the mean field model
for the yielding transition.

Figure 4a shows the unified yielding state diagram for soft colloids obtained
using reduced variables. For each sample, we tune $\bar{\alpha}$ and the
variance $\sigma^{2}_{\alpha}$ of the coupling constants distribution, the
spontaneous relaxation rate
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}1/\tau_{0}}}$
and the exponent $n$ in order to reproduce the strain dependence of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{s,f}}}$
and $\chi$, as exemplified in Fig. 2d-f. Using these fitting parameters, we
re-express the experimental variables in terms of the reduced variables of
Fig. 4a. In this representation, all samples fall on nearly the same solid and
fluid branches, characterized respectively by a single, compressed exponential
slow mode (small
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}}}$,
blue solid symbols in Fig. 4a) and a single, fast stretched exponential mode
(large
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}}}$,
red open symbols). This collapse is remarkable, given the diversity of the
microscopic structure of the investigated samples, which in turn results in
marked differences in the sensitivity to shear, compare e.g. the steep growth
of
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{s}}}$
with applied strain for the emulsions to the gentler increase for microgels
and nanoparticles (Figs. 2d-f). Furthermore, by analyzing data at various $q$
vectors for E samples, we find that the collapse is robust with respect to the
choice of the probed length scale, see Supplementary. At intermediate
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}}}$,
within the region inaccessible to the mean field model, a fast mode and a slow
mode coexist (gray half-filled symbols in Fig. 4a), as predicted by the model
with disorder; the abscissa used for these points is the weighted average of
the fast and slow relaxation rates (see Methods).

In the coexisting region, samples’ properties vary markedly and systematically
with $K_{r}$, which suggests classifying all systems according to this
parameter. Since $K_{r}<1$ corresponds to glassy samples and $K_{r}>1$ to
equilibrated fluids, we quantify ‘glassiness’ of samples with $K_{r}<1$ by
$g=1-K_{r}$, which increases as the iso-$K_{r}$ curves move downward (see
arrow in Fig. 4a) away from the critical curve ($g=0$, $K_{r}=1$). For samples
of the same kind, the trend in $g$ (see Fig. 4b and Table SI2 in
Supplementary) agrees with the behavior intuitively associated with a lesser
or greater glassiness. For emulsions and microgels, we find that in general
the higher $\varphi$ the more glassy the sample. Consistent with the notion
that with age systems evolve towards deeper states in the glassy phase, we
find that $g$ increases with age for sample N41%. Finally, one expects $g$ to
increase with $\omega$, because glassy samples fall increasingly out of
equilibrium as the time scale of the driving becomes shorter. This is indeed
what is seen in our experiments, compare samples M40s and M2s. Beyond these
comparisons, the notion of glassiness paves the way for a quantitative
comparison between samples of different nature (e.g. emulsions vs microgels)
or probed according to different protocols. Keeping in mind that each
iso-$K_{r}$ curve –and thus each $g$ value– is characterized by a different
disorder parameter of the coupling constants
$d=\sigma_{\alpha}^{2}/\bar{\alpha}^{2}$, we find a remarkable negative
correlation between the glassiness and $d$, which, as shown in Fig. 4b, is
characterized by a simple master curve.

Figure 4:  Unified state diagram for the yielding transition. a): Yielding
state diagram for all samples, using reduced variables (see text). Blue, solid
(resp., red, open) symbols: samples with a single slow (resp., fast)
relaxation mode. Gray half-filled symbols: coexistence of the slow and fast
modes. The symbol shape is the same as in b), where the samples are identified
by labels, see Methods. The continuous, dotted and dash-dotted lines indicate
the solid, fluid and coexistence branches of the iso-$K$ lines obtained from
simulations of the model with disorder (same color codes as in b)). Area
shaded in light blue below the black line: region inaccessible to the mean
field model, corresponding to the coexistence of fast and slow relaxation
modes in experiments and simulations of the model with disorder. b): Disorder
$d=\sigma_{\alpha}^{2}/\bar{\alpha}^{2}$ of the coupling constants
$\alpha_{i,j}$ (see Eq. 2) as a function of ‘glassiness’, defined as the
reduced distance form the critical point, $g=1-K_{r}$ (see also the arrow in
a)).

The correlation between glassiness and disorder has also deep implications on
the nature of the yielding transition. We find that the most glassy samples,
for which $d$ approaches the $d=0$ mean field limit, exhibit features typical
of a first-order transition, as predicted by the mean field model (Fig. 5a-c).
Figure 5a shows $w$, the relative width of the $\gamma_{0}$ range where fast
and slow modes coexists in experiments, demonstrating a dramatic increase of
the transition sharpness as $d$ decreases. Figures 5b-c demonstrate
hysteresis, another distinctive feature of first-order transitions (data from
simulations). As seen in Fig. 5c, hysteresis is largest for the smallest
disorder and vanishes when departing from the mean field limit.

Conversely, systems close to the critical $K$ and hence with large $d$ (small
$g$) exhibit features usually associated with second order transitions,
illustrated in Fig. 5d-f. Figure 5d shows for sample M2s that temporal
fluctuations of the relaxation rate are strongly enhanced at the transition.
Figures 5e-f demonstrate sluggishness, another feature of second-order
transitions: in both experiments (Fig. 5e) and simulations (Fig. 5f) on
systems with large $d$ (small $g$), the time to attain a stationary state
dramatically increases around the transition.

Figure 5:  Distance from a critical point governs the nature of the yielding
transition. Top row: features typical of a first order transition are enhanced
as glassiness increases ($d$ decreases). a): Relative width $w$ of the
transition in experiments probing the microscopic dynamics across yielding
(see Methods for details). b): Hysteresis at yielding, data from simulations
with $d=0$ (dashed lines, blue shades) and $d=0.02$ (solid lines, red shades).
c): Area $H$ of the hysteresis loop in simulations with disorder, normalized
by that in the $d=0$ limit. Bottom row: features typical of a second order
transition emerge close to the critical point (small $g$, large $d$). d):
Fluctuations of the dynamics, quantified by the temporal variance
$\sigma^{2}_{\tau_{r}}$ of the relaxation time of the intensity correlation
functions, are strongly enhanced at the transition. Data for sample M2s. e):
Normalized time $t_{eq}\omega$ to attain stationary dynamics for sample N41%.
The lines are a guide to the eye. f): Number of iterations required for
numerical solutions of the model with disorder to converge. Data are
normalized by their corresponding value at $d=0$. Black symbols and line:
$d=0.1$; blue line: $d=10^{-4}$. In d)-f), the gray shaded area indicates the
range of the dynamical transition (see Methods for details). Panels b, c, f:
same mean field parameters as in Figs. 2b,e.

The unified yielding state diagram established here demonstrates the universal
nature of the solid-fluid transition in soft glassy systems under oscillatory
shear. Our model differs from previous approaches by introducing a direct link
between the macroscopic drive and the microscopic dynamics, and by including
in the latter the contribution of ultra-slow, spontaneous relaxations; as a
result, slow and fast relaxation modes may coexist. These two modes might be
associated with the formation of shear bands organized in the shear gradient
direction [42, 43]. Our light scattering experiments cannot test directly this
hypothesis; however, our numerical data, microscopy experiments on emulsions,
and photon correlation imaging measurements on sample N45% are inconsistent
with shear banding (Supplementary Sec. IV), rather pointing to the formation
of fast-relaxing domains similar to those of Fig. 3c. This suggests that the
coexistence of slow- and fast-relaxing regions is a generic feature, not
necessarily implying the structuration in shear bands. A key result of our
work is the emergence of different features in the transition depending on the
distance from the critical point. This indicates a promising way to reconcile
apparently conflicting reports in previous studies of yielding.

Finally, the model proposed here does not depend on the details of the
interaction between the microscopic constituents of the system: we thus expect
it to provide a general framework for the yielding transition, possibly
including in systems with attractive interactions [44, 45].

## Methods

### Samples

Table 2 summarizes the main sample and experimental parameters. PNIPAM
microgels (M) were synthesized by emulsion polymerization as described in
[46], and were suspended in a 2 mM $\mathrm{NaN_{3}}$ aqueous solution to
prevent bacterial growth. The microgel radius at $T=23~{}^{\circ}\mathrm{C}$
is 294 nm, as measured by dynamic light scattering (DLS) in a diluted
suspension. Sample preparation and characterization, including the
determination of the effective volume fraction $\varphi$, are described in
[47]. Note that the effective volume fraction is larger than one, due to the
particle softness. Charged nanoparticle systems (N) were prepared by
concentrating an acqueous suspension of silica particles (Ludox TM50, from
Sigma Aldrich), as described in [47]. The particles have a hydrodynamic radius
of 23 nm, as measured by DLS in the dilute limit. To improve the scattering
signal, the samples were seeded with 200 nm-diameter polystyrene particles at
extremely low volume fraction, $\phi_{PS}<10^{-6}$. We check that this seeding
has no measurable impact on the rheological properties microscopic dynamics of
the samples. Concentrated emulsions (E) were prepared by dispersing
polydimethylsiloxane droplets into a water/glycerol matrix, as described in
[4]. The resulting droplets have an average size of $2.4~{}\mu\mathrm{m}$ and
20% polydispersity. All samples were initialized by applying a preshear, see
Sec. IIa of the Supplementary for details. Rheology and microscopic dynamics
started immediately after applying the preshear for all samples, except for
sample N41%Aged, which was left at rest for 12h prior to measurements.

Sample ID | $\varphi$ |  $T$ (s) |  $\omega$ [rad/s] |  $q$ $[~{}\mu\mathrm{m}^{-1}]$ | setup |  $\mathbf{q}$ orientation
---|---|---|---|---|---|---
M2s | 1.5 | 2 | 3.14 | 5 | PCI | vorticity
M40s | 1.5 | 40 | 0.157 | 0.1-5 | SALS | vorticity
N41% | 0.41 | 2 | 3.14 | 30 | PCI | shear gradient
N41%Aged | 0.41 | 2 | 3.14 | 30 | PCI | shear gradient
N45% | 0.45 | 2 | 3.14 | 30 | PCI | shear gradient
E65% | 0.65 | 1 | 6.28 | 1-20 | ff-DDM | vorticity
E70% | 0.70 | 1 | 6.28 | 1-20 | ff-DDM | vorticity
E74% | 0.74 | 1 | 6.28 | 1-20 | ff-DDM | vorticity
E88% | 0.88 | 1 | 6.28 | 1-20 | ff-DDM | vorticity
Table 2: Main features of the samples used in this study. M, N, and E refer to
microgels, nanoparticles and emulsions, respectively. $\varphi$ designates the
volume fraction for N and E, and the effective volume fraction for M, as
defined in [47]. $T$ is the period of the oscillatory strain. PCI, SALS, and
ff-DDM are photon correlation imaging, small angle light scattering and far
field differential dynamic microscopy, respectively. The last column shows the
orientation of the largest component of $\mathbf{q}$ with respect to the shear
field, see Experimental setups for details.

### Experimental setups

With the only exception of E samples, all experiments are performed with a
home-made shear cell equipped with sliding parallel plates [48], sketched in
Supplementary Fig. SI2. One plate is driven by a piezoelectric strain actuator
(P602, from Physik Instrumente) and a force sensor (LC601, from Omega
Engineering) measures the force applied by the actuator, enabling strain-
controlled rheology experiments. The sample has a cross-sectional area of
about $4~{}\mathrm{cm}^{2}$ and a thickness between $300$ and
$700~{}~{}\mu\mathrm{m}$. For most experiments, shear rheology was coupled to
a spatially-resolved Photon Correlation Imaging (PCI) apparatus [49],
collecting light scattered in a direction orthogonal to the shear and forming
a scattering angle $\theta$ with the incoming beam. For samples N, we choose
$\theta=141^{\circ}$, yielding a scattering vector $q=4\pi
n_{r}\lambda^{-1}\sin(\theta/2)=30~{}\mu\mathrm{m}^{-1}$, where $n_{r}=1.33$
is the solvent refractive index and $\lambda=532.5$ nm the wavelength of laser
light. In this case, $q$ is predominantly oriented along the shear gradient
direction, with a minor component $q_{x}=3~{}\mu\mathrm{m}^{-1}$ along the
vorticity direction. For M2s, we choose $\theta=20^{\circ}$ such that
$q=5~{}\mu\mathrm{m}^{-1}$ is predominantly oriented along the vorticity
direction, with a minor component $q_{z}=0.9~{}\mu\mathrm{m}^{-1}$ along the
shear gradient. For M40s, experiments are performed using a different
scattering geometry: a far-field small angle light scattering apparatus (SALS)
enabling multiple scattering vectors to be probed simultaneously, oriented
along both the velocity and the vorticity direction, with
$0.1~{}~{}\mu\mathrm{m}^{-1}\leq q\leq 5~{}~{}\mu\mathrm{m}^{-1}$ . Data shown
in the main text correspond to $q=4.8~{}~{}\mu\mathrm{m}^{-1}$, oriented along
the vorticity direction. See Supplementary Figs SI1-SI3 for the schemes of the
setups.

For samples E, a different home-made shear cell with parallel, counter-
translating plates displaced by a piezoelectric actuator is coupled to an
inverted microscope with differential interference contrast (DIC) optics [4] ,
see Supplementary Fig. SI3. The acquired imaged region has a depth of field of
$0.5~{}~{}\mu\mathrm{m}$, much smaller than the droplet size, and the imaged
plane corresponds to the stagnation plane of the shear deformation. We analyze
microscopy videos using far-field Differential Dynamic Microscopy (ff-DDM)
[50], which yields an intensity correlation function equivalent to DLS, with
scattering vectors $1~{}~{}\mu\mathrm{m}^{-1}\leq q\leq
20~{}~{}\mu\mathrm{m}^{-1}$ in the velocity-vorticity plane. Data presented in
the main text correspond to $q=9~{}~{}\mu\mathrm{m}^{-1}$ along the vorticity
direction. Data for more scattering vectors are included in Supplementary Fig.
SI6-SI9.

### Characterization of the microscopic dynamics

We quantify the microscopic dynamics via the two-time intensity correlation
$g_{2}(t,\tau)-1=\beta\langle I_{p}(t)I_{p}(t+\tau)\rangle/[\langle
I_{p}(t)\rangle\langle I_{p}(t+\tau)\rangle]$, with $\beta$ a constant such
that $g_{2}(t,\tau)-1\rightarrow 1$ for $\tau\rightarrow 0$. $I_{p}$ is the
scattered intensity collected by the $p-th$ pixel of the detector (for the PCI
and SALS setups), or a component of the Fourier transform of the microscope
images for sample E. $\langle\dots\rangle$ indicates the average over a set of
pixels corresponding nearly to the same scattering vector or Fourier
component.

In the stationary regime, we average $g_{2}-1$ over time $t$ to reduce noise
before fitting Eq. 1 to the data. In Fig. 4a, the abscissa of the state points
belonging to the coexistence region is calculated as the weighted average of
the slow and fast relaxation rates, normalized by the relaxation rate at the
critical point:
${{\color[rgb]{0,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{0,0,0}\pgfsys@color@gray@stroke{0}\pgfsys@color@gray@fill{0}\Gamma_{r}=[\chi\Gamma_{s}+(\chi-1)\Gamma_{f}]/\Gamma_{c}}}$.
The slow and fast relaxation rates are obtained by fitting Eq. 1 to the
experimental $g_{2}-1$.

To quantify the temporal fluctuations of the dynamics, we inspect the
$t$-dependence of the two-time correlation functions, with no $t$ averaging
performed on $g_{2}(t,\tau)-1$. For each $t$, a relaxation time $\tau_{r}(t)$
is obtained from $\tau_{r}=\int_{0}^{\infty}[g_{2}(t,\tau)-1]\mathrm{d}\tau$,
a procedure more robust than fitting to Eq. 1 when dealing with correlation
functions that are not averaged over $t$ and are thus more noisy. In Fig. 5d
we show $\sigma^{2}_{\tau_{r}}$, the temporal variance of $\tau_{r}(t)$.

We locate the range of the dynamic transition by determining the two strain
amplitudes $\gamma_{0,s}$ and $\gamma_{0,e}$ such that $0.05\leq\chi\leq 0.95$
for $\gamma_{0,s}\leq\gamma_{0}\leq\gamma_{0,e}$. The gray shaded regions in
Figs. 5d-f highlight the strain range
$\gamma_{0,s}\leq\gamma_{0}\leq\gamma_{0,e}$. The normalized width of the
transition shown in Fig. 5a is defined as
$w=(\gamma_{0,e}-\gamma_{0,s})/(\gamma_{0,e}+\gamma_{0,s})$.

### Numerical solution of the model with disorder

To study the effect of disorder on yielding, we implement our model (Eq. 2) on
a $D$-dimensional cubic lattice with periodic boundary conditions. Each site,
$i$, is assigned a local relaxation rate, $\Gamma_{i}$, and each pair of
neighbor sites is attributed a coupling constant, $\alpha_{ij}$, randomly
drawn from a probability distribution, $P(\alpha)$. For a given strain
amplitude $\gamma_{0}$ and starting from an initial configuration of local
rates $\Gamma_{i}^{(init)}$, we seek a configuration of local rates
$\Gamma_{i}^{(sol)}$ that satisfies Eq. 2 for all sites. In our
implementation, $\Gamma_{i}^{(sol)}$ is approached iteratively: at each step,
a set of target site rates $\Gamma_{i}^{(t)}$ is computed through Eq. 2 using
the current site rates $\Gamma_{i}^{(c)}$. $\Gamma_{i}^{(t)}$ then replaces
$\Gamma_{i}^{(c)}$ for the following iteration, yielding a new set of target
site rates. The convergence criterion is expressed in terms of a loss function
$\mathcal{L}=\sum_{i}\left[\left(\log\Gamma_{i}^{(c)}-\log\Gamma_{i}^{(t)}\right)/\log\Gamma_{i}^{(t)}\right]^{2}$,
which tends to 0 as $\Gamma_{i}^{(c)}$ approaches $\Gamma_{i}^{(sol)}$.

Results shown in this paper are obtained for $D=2$ and a Log-Normal
distribution of the coupling constants:
$P(\alpha;\mu,\sigma)=(\alpha\sigma\sqrt{2\pi})^{-1}\exp[-(\ln\alpha-\mu)^{2}/(2\sigma^{2})]$,
where the average $\bar{\alpha}$ and variance $\sigma_{\alpha}^{2}$ of the
coupling constants are related to the parameters $\mu,\sigma$ of the Log-
Normal distribution by $\bar{\alpha}=\exp(\mu+\sigma^{2}/2)$ and
$\sigma_{\alpha}^{2}=\bar{\alpha}^{2}(\exp\sigma^{2}-1)$. We quantify disorder
by the dimensionless parameter $d=\sigma_{\alpha}^{2}/\bar{\alpha}^{2}$.
Representative results for different choices of $P(\alpha)$ are shown in
Supplementary Fig. SI13, and exhibit no qualitative differences.

To mimic the effect of preshear in experiments, the iterative solution of the
numerical model is typically initiated with a uniform distribution of local
rates $\Gamma_{i}^{(init)}=\Gamma^{(0)}$. The effect of hysteresis shown in
Figs. 5b,c is studied by simulating an amplitude sweep experiment:
$\gamma_{0}$ is first increased from low to high amplitudes and then decreased
from high to low amplitudes, every time initiating the iterative calculation
from the model solution for the previous amplitude.

## Data availability

The data that support the plots within this paper and other findings of this
study are available from the corresponding authors upon request.

## Code availability

The code that support the plots within this paper and other findings of this
study are available from the corresponding authors upon request.

## References

  * Bonn  _et al._ [2017] D. Bonn, M. M. Denn, L. Berthier, T. Divoux, and S. Manneville, Yield Stress Materials in Soft Condensed Matter, Rev. Mod. Phys. 89, 035005 (2017).
  * Koumakis  _et al._ [2013] N. Koumakis, J. F. Brady, and G. Petekidis, Complex Oscillatory Yielding of Model Hard-Sphere Glasses, Physical Review Letters 110, 178301 (2013).
  * Mason  _et al._ [1996] T. G. Mason, J. Bibette, and D. A. Weitz, Yielding and flow of monodisperse emulsions, Journal of Colloid and Interface Science 179, 439 (1996).
  * Knowlton  _et al._ [2014] E. D. Knowlton, D. J. Pine, and L. Cipelletti, A microscopic view of the yielding transition in concentrated emulsions, Soft Matter 10, 6931 (2014).
  * Rogers  _et al._ [2018] M. C. Rogers, K. Chen, M. J. Pagenkopp, T. G. Mason, S. Narayanan, J. L. Harden, and R. L. Leheny, Microscopic signatures of yielding in concentrated nanoemulsions under large-amplitude oscillatory shear, Physical Review Materials 2, 095601 (2018).
  * Rogers  _et al._ [2011] S. A. Rogers, B. M. Erwin, D. Vlassopoulos, and M. Cloitre, Oscillatory yielding of a colloidal star glass, Journal of Rheology 55, 733 (2011).
  * Ketz  _et al._ [1988] R. J. Ketz, R. K. Prud’homme, and W. W. Graessley, Rheology of concentrated microgel solutions, Rheol Acta 27, 531 (1988).
  * Sollich  _et al._ [1997] P. Sollich, F. Lequeux, P. Hébraud, and M. E. Cates, Rheology of soft glassy materials, Physical Review Letters 78, 2020 (1997).
  * Seth  _et al._ [2011] J. R. Seth, L. Mohan, C. Locatelli-Champagne, M. Cloitre, and R. T. Bonnecaze, A micromechanical model to predict the flow of soft particle glasses, Nature Materials 10, 838 (2011).
  * Donley  _et al._ [2020] G. J. Donley, P. K. Singh, A. Shetty, and S. A. Rogers, Elucidating the G′′ overshoot in soft materials with a yield transition via a time-resolved experimental strain decomposition, Proc Natl Acad Sci USA 117, 21945 (2020).
  * Brader  _et al._ [2010] J. M. Brader, M. Siebenbürger, M. Ballauff, K. Reinheimer, M. Wilhelm, S. J. Frey, F. Weysser, and M. Fuchs, Nonlinear response of dense colloidal suspensions under oscillatory shear: Mode-coupling theory and fourier transform rheology experiments, Physical Review E 82, 061401 (2010).
  * Voigtmann [2014] T. Voigtmann, Nonlinear glassy rheology, Current Opinion in Colloid & Interface Science 19, 549–560 (2014).
  * Picard  _et al._ [2002] G. Picard, A. Ajdari, L. Bocquet, and F. Lequeux, Simple model for heterogeneous flows of yield stress fluids, Physical Review E 66, 051501 (2002).
  * Benzi  _et al._ [2019] R. Benzi, T. Divoux, C. Barentin, S. Manneville, M. Sbragaglia, and F. Toschi, Unified Theoretical and Experimental View on Transient Shear Banding, Phys. Rev. Lett. 123, 248001 (2019).
  * Liu  _et al._ [2018] C. Liu, K. Martens, and J.-L. Barrat, Mean-Field Scenario for the Athermal Creep Dynamics of Yield-Stress Fluids, Phys. Rev. Lett. 120, 028004 (2018).
  * Sainudiin  _et al._ [2015] R. Sainudiin, M. Moyers-Gonzalez, and T. Burghelea, A microscopic Gibbs field model for the macroscopic yielding behaviour of a viscoplastic fluid, Soft Matter 11, 5531 (2015).
  * Keim and Arratia [2013] N. C. Keim and P. E. Arratia, Yielding and microstructure in a 2D jammed material under shear deformation, Soft Matter 9, 6222 (2013).
  * Fiocco  _et al._ [2013] D. Fiocco, G. Foffi, and S. Sastry, Oscillatory athermal quasistatic deformation of a model glass, Physical Review E 88, 020301 (2013).
  * Hima Nagamanasa  _et al._ [2014] K. Hima Nagamanasa, S. Gokhale, A. K. Sood, and R. Ganapathy, Experimental signatures of a nonequilibrium phase transition governing the yielding of a soft glass, Physical Review E 89, 10.1103/PhysRevE.89.062308 (2014).
  * Jeanneret and Bartolo [2014] R. Jeanneret and D. Bartolo, Geometrically protected reversibility in hydrodynamic Loschmidt-echo experiments, Nat Commun 5, 3474 (2014).
  * Kawasaki and Berthier [2016] T. Kawasaki and L. Berthier, Macroscopic yielding in jammed solids is accompanied by a nonequilibrium first-order transition in particle trajectories, Physical Review E 94, 022615 (2016).
  * Leishangthem  _et al._ [2017] P. Leishangthem, A. D. S. Parmar, and S. Sastry, The yielding transition in amorphous solids under oscillatory shear deformation, Nat Commun 8, 14653 (2017).
  * Edera  _et al._ [2021] P. Edera, M. Brizioli, G. Zanchetta, G. Petekidis, F. Giavazzi, and R. Cerbino, Deformation profiles and microscopic dynamics of complex fluids during oscillatory shear experiments, Soft Matter 17, 8553–8566 (2021).
  * Lerouge and Berret [2009] S. Lerouge and J.-F. Berret, Shear-Induced Transitions and Instabilities in Surfactant Wormlike Micelles, in  _Polymer Characterization_, Vol. 230, edited by K. Dusek and J.-F. Joanny (Springer Berlin Heidelberg, Berlin, Heidelberg, 2009) pp. 1–71.
  * Cipelletti  _et al._ [2003] L. Cipelletti, L. Ramos, S. Manley, E. Pitard, D. A. Weitz, E. E. Pashkovski, and M. Johansson, Universal non-diffusive slow dynamics in aging soft matter, Faraday Discuss. 123, 237 (2003).
  * Madsen  _et al._ [2010] A. Madsen, R. L. Leheny, H. Guo, M. Sprung, and O. Czakkel, Beyond simple exponential correlation functions and equilibrium dynamics in x-ray photon correlation spectroscopy, New Journal of Physics 12, 055001 (2010).
  * Denisov  _et al._ [2015] D. V. Denisov, M. T. Dang, B. Struth, A. Zaccone, G. H. Wegdam, and P. Schall, Sharp symmetry-change marks the mechanical failure transition of glasses, Scientific Reports 5, 14359 (2015).
  * Divoux  _et al._ [2013] T. Divoux, V. Grenard, and S. Manneville, Rheological Hysteresis in Soft Glassy Materials, Physical Review Letters 110, 018304 (2013).
  * Bocquet  _et al._ [2009] L. Bocquet, A. Colin, and A. Ajdari, Kinetic Theory of Plastic Flow in Soft Glassy Materials, Physical Review Letters 103, 036001 (2009).
  * Nordstrom  _et al._ [2011] K. N. Nordstrom, J. P. Gollub, and D. J. Durian, Dynamical heterogeneity in soft-particle suspensions under shear, Physical Review E 84, 021403 (2011).
  * Hebraud  _et al._ [1997] P. Hebraud, F. Lequeux, J. P. Munch, and D. J. Pine, Yielding and Rearrangements in Disordered Emulsions, Phys. Rev. Lett. 78, 4657 (1997).
  * van Megen  _et al._ [1998] W. van Megen, T. C. Mortensen, S. R. Williams, and J. Muller, Measurement of the self-intermediate scattering function of suspensions of hard spherical particles near the glass transition, Phys. Rev. E 58, 6073 (1998).
  * Weeks  _et al._ [2000] E. R. Weeks, J. C. Crocker, A. C. Levitt, A. Schofield, and D. A. Weitz, Three-dimensional direct imaging of structural relaxation near the colloidal glass transition, Science 287, 627 (2000).
  * Derec  _et al._ [2001] C. Derec, A. Ajdari, and F. Lequeux, Rheology and aging: A simple approach, The European Physical Journal E: Soft Matter and Biological Physics 4, 355 (2001).
  * Miyazaki  _et al._ [2006] K. Miyazaki, H. M. Wyss, D. A. Weitz, and D. R. Reichman, Nonlinear viscoelasticity of metastable complex fluids, Europhysics Letters (EPL) 75, 915 (2006).
  * Hess and Aksel [2011] A. Hess and N. Aksel, Yielding and structural relaxation in soft materials: Evaluation of strain-rate frequency superposition data by the stress decomposition method, Physical Review E 84, 051502 (2011).
  * Biroli and Garrahan [2013] G. Biroli and J. P. Garrahan, Perspective: The glass transition, The Journal of Chemical Physics 138, 12A301 (2013).
  * Zausch  _et al._ [2008] J. Zausch, J. Horbach, M. Laurati, S. U. Egelhaaf, J. M. Brader, T. Voigtmann, and M. Fuchs, From equilibrium to steady state: the transient dynamics of colloidal liquids under shear, Journal of Physics: Condensed Matter 20, 404210 (2008).
  * [39] S. Aime et al., Microscopic yielding of glassy materials under oscillatory shear.
  * Berker [1993] A. N. Berker, Critical behavior induced by quenched disorder, Physica A: Statistical Mechanics and its Applications 194, 72–76 (1993).
  * Bellafard  _et al._ [2015] A. Bellafard, S. Chakravarty, M. Troyer, and H. G. Katzgraber, The effect of quenched bond disorder on first-order phase transitions, Annals of Physics 357, 66–78 (2015).
  * Divoux  _et al._ [2016] T. Divoux, M. A. Fardin, S. Manneville, and S. Lerouge, Shear banding of complex fluids, Annual Review of Fluid Mechanics 48, 81–103 (2016).
  * Radhakrishnan and Fielding [2016] R. Radhakrishnan and S. M. Fielding, Shear Banding of Soft Glassy Materials in Large Amplitude Oscillatory Shear, Physical Review Letters 117, 188001 (2016).
  * Pham  _et al._ [2002] K. N. Pham, A. M. Puertas, J. Bergenholtz, S. U. Egelhaaf, A. Moussaid, P. N. Pusey, A. B. Schofield, M. E. Cates, M. Fuchs, and W. C. K. Poon, Multiple glassy states in a simple model system, Science 296, 104 (2002).
  * Gibaud  _et al._ [2010] T. Gibaud, D. Frelat, and S. Manneville, Heterogeneous yielding dynamics in a colloidal gel, Soft Matter 6, 3482 (2010).
  * Truzzolillo  _et al._ [2018] D. Truzzolillo, S. Sennato, S. Sarti, S. Casciardi, C. Bazzoni, and F. Bordi, Overcharging and reentrant condensation of thermoresponsive ionic microgels, Soft Matter 14, 4110 (2018).
  * Philippe  _et al._ [2018] A.-M. Philippe, D. Truzzolillo, J. Galvan-Myoshi, P. Dieudonné-George, V. Trappe, L. Berthier, and L. Cipelletti, Glass transition of soft colloids, Physical Review E 97, 040601(R) (2018).
  * Aime  _et al._ [2016] S. Aime, L. Ramos, J. M. Fromental, G. Prévot, R. Jelinek, and L. Cipelletti, A stress-controlled shear cell for small-angle light scattering and microscopy, Review of Scientific Instruments 87, 123907 (2016).
  * Cipelletti  _et al._ [2016] L. Cipelletti, V. Trappe, and D. J. Pine, Scattering Techniques, in  _Fluids, Colloids and Soft Materials_ , edited by A. Fernandez-Nieves and A. Puertas (John Wiley & Sons, Inc., 2016) pp. 131–148.
  * Aime and Cipelletti [2019] S. Aime and L. Cipelletti, Probing shear-induced rearrangements in Fourier space. II. Differential dynamic microscopy, Soft Matter 15, 213 (2019).

###### Acknowledgements.

We thank E. D. Knowlton for help with the experiments on emulsions and L.
Berthier for illuminating discussions. This work was funded by the French
CNES, ANR (grants No. ANR-14-CE32-0005, FAPRES, and ANR-20-CE06-0028,
MultiNet), and by the EU (Marie Sklodowska-Curie ITN Supolen Grant 607937). LC
acknowledges support from the Institut Universitaire de France.

## Author contributions

SA, LR, DJP, and LC designed experiments. SA performed experiments and
numerical simulations. SA and DT conceived the model. All authors analyzed the
results, discussed and improved the model, and contributed to writing the
paper.

## Competing interests

The authors declare no competing interests.

## Additional information

Supplementary information The online version contains supplementary material
available at [url to be inserted].

Correspondence and requests for materials should be addressed to SA or LC.