Nonlinear rheological behavior of glass-forming colloidal suspensions under oscillatory shear: Experiment and relation to mode coupling theory predictions

Nano- or micrometer-sized solid particles dispersed in a liquid have extensive technological applications, ranging from cosmetics and food science to building materials. Most such colloidal suspensions use water as the continuous phase. Examples include toothpaste, ink, paints, blood, or clay. These examples demonstrate the relevance of colloidal suspensions to a broad spectrum of scientific research [1]. Understanding the flow behavior of colloidal suspensions is a critical part of adjusting a particular colloidal system for specific applications, for instance, when mechanical behavior must be adjusted or enhanced [2]. When the colloid volume fraction becomes appreciable, the viscosity of the suspension increases and it exhibits non-Newtonian mechanical behavior as one of its characteristics [1,2].

If the particles in a colloidal system are close to monodisperse and spherical, they can either form a crystalline state or become trapped in a metastable glassy state. Which state is formed depends on how fast the critical volume fraction is reached. If the critical volume fraction is reached quickly, the crystallization process is kinetically suppressed, leading to a glassy state. Shear induced crystallization can occur for colloids in the metastable glassy state [3]. To investigate a pure glassy state, the crystallization has to be prevented by introducing a certain dispersity of the radius of the colloids (⁠ $σ r e l>12%$ [4]). In general, the glassy state of colloids is characterized by trapping the particles in a cage formed by the surrounding particles in which the movement of the particles is cooperatively constrained [5]. This constraint of particle movement causes the colloidal suspension to exhibit viscoelastic behavior in the glassy state as well as a yield stress. These topics have been extensively investigated [6–14].

In addition, colloidal suspensions exhibit significant nonlinear mechanical behavior when exposed to large deformation, as shown in several publications [13,15–20]. Understanding this nonlinear behavior is crucial since materials frequently undergo large deformations during industrial processing and real-life applications.

At large strain amplitudes, the material's response functions, namely, the storage modulus $G ′$ and the loss modulus $G ′ ′$⁠, begin to display a dependence on the applied strain amplitude. For most materials, this is accompanied by an anharmonic stress response of the material, meaning that the stress response cannot be described by a single sine wave. Just a few materials show only a dependence of the moduli on the applied strain amplitude while still having a pure sinusoidal stress response [21–26].

There are several ways to analyze the anharmonic stress response. The first approach was plotting $σ(t)$ vs $γ(t)$ or vs $γ ˙(t)$⁠, obtaining the so-called Lissajous–Bowditch figures, while the first quantitative approach to analyze the anharmonicity of a material was Fourier transform (FT) analysis [19,27–41]. With this approach, the higher harmonics of the stress signal are obtained. Another method that has been developed is a stress decomposition into elastic and viscous components [42,43], which gives some physical meaning to the anharmonic behavior, or to analyze $σ ′(t)$ and $σ ′ ′(t)$ as a series of Chebyshev polynomials [44–47], which also makes a distinction between the elastic and the viscous contributions to the anharmonicity. All those methods are based on the full cycle of the oscillatory response. As an alternative to those, the sequence of physical process (SPP) analysis has been suggested [48–51], where the transient stress response of the material is analyzed having the advantage over the other two methods of enabling microstructural interpretations throughout the oscillation cycle [52,53]. However, FT rheology and the Chebyshev method have the advantage that they are easy to calculate [54].

Recently, recovery rheology was developed [55–57], which can be used to understand the full amplitude sweep by splitting the strain in its recoverable and unrecoverable components. Utilizing recovery rheology [57] by combining strain-controlled oscillatory shear with stress-controlled recovery tests, it can be distinguished between recoverable and unrecoverable contributions to the total energy dissipation. Furthermore, traditional material functions can be redefined using recovery rheology unifying the interpretations of MAOS response [54].

The focus of this work is to answer a question raised by a previous comparison of FT rheological results to mode coupling theory (MCT). There, the strong growth of the third harmonic response when approaching the glass transition was observed in LAOS strain measurements, while theory even predicts its divergence [58]. The behavior in the glass state was not addressed. Therefore, FT rheology is chosen to gain a quantitative analysis of the nonlinear response approaching the colloidal glass transition and in the colloidal glass which enables the comparison to the MCT calculations.

Using the FT approach, the nonlinear stress signal is expressed as a Fourier series as the stress response can be interpreted as a sum of sines with the fundamental frequency of the applied oscillatory strain and the odd harmonics of the fundamental frequency with phase shifts. Because of symmetry reasons and, as described in detail in Wilhelm [30], Cziep et al. [59], and Hirschberg et al. [60], the even harmonics are not present as long as flow instabilities such as slip do not occur [20]. Note that the FT-rheology spectra are initially recorded when analyzing data via Chebyshev polynomials. Furthermore, in the analysis via the SPP approach, the FT can be applied as an efficient comb filter in frequency space.

The relative intensity of the magnitude of the third harmonic with respect to the fundamental, $I 3 / 1= I 3/ I 1$⁠, obtained from the FT is utilized as a quantification of the anharmonicity of the investigated material [41] and can be used as a quantitative measure for the transition from the linear region to the asymptotic nonlinear region of the material, the so-called medium amplitude oscillatory shear (MAOS) region [61–63]. Another measure of the onset of nonlinearities, which is often used, is the asymptotic deviation of $G ′$ from its plateau value. The advantage of using $I 3 / 1$ is the possibility to extrapolate it to a specific threshold, e.g., $I 3 I 1<0.005$⁠, to gain a characteristic strain amplitude for the onset of nonlinear behavior. Furthermore, the anharmonicity $I 3$ can be quantified at the third harmonic frequency $3 ω 1$⁠. At this frequency, the stress response vanishes in linear response, and $I 3$ gives the dominant anharmonic response. It decreases for decreasing strain amplitude according to scaling laws.

However, only a few oscillatory nonlinear rheological measurements have been performed on glass-forming colloidal suspensions, and most are limited to a few frequencies and focus on the dependence of the anharmonicity on the strain amplitude [13,64–66]. Only a few publications address the frequency dependence of the anharmonic behavior [15,16,58,67]. In Seyboldt et al. [58], the frequency dependence of the anharmonicity is addressed in the supercooled state of the sample close to the glass transition volume fraction and compared to MCT. The MCT predictions show that the anharmonicity of the colloidal suspension strongly grows in an intermediate frequency window when approaching the glass transition volume fraction. Upon increasing the volume fraction, this maximum increases and shifts to lower frequencies, with a predicted power-law divergence for low frequencies at the glass transition. This behavior raises the question of what is happening within the glass. Preliminary data in the glass already exist [16], but extensive experimental data, especially in the low-frequency regime, are still missing.

To fill this gap, we scaled up an existing synthesis of colloidal suspensions of polystyrene-poly-N-isopropylacrylamide core-shell particles [68] to be able to extensively study the frequency dependence of the anharmonicity of colloidal suspensions in the fluidlike and the glassy regimes. In contrast to previous investigations, which used monomodal samples with a high dispersity in radius to suppress the particle crystallization process, we study a well characterized binary mixture of spheres, whose individual radius distributions are rather narrow for gaining better comparability to MCT calculations as they are mostly based on monodisperse samples and can be extended to binary mixtures [69].

In this article, we present oscillatory nonlinear measurements at multiple frequencies on this model system below and above the glass transition volume fraction. The core-shell particles possess an intermediate interparticle (soft) potential, in contrast to the starlike micelles (ultrasoft) and the sterically stabilized polymethyl methacrylate (PMMA) hard spheres (hard) used as the two extrema studied in Poulos et al. [16]. We compare the strain amplitude dependence of $I 3 / 1$ in the glass and the liquidlike region close to the glass transition of our soft system with the results of Poulos et al. [16] on ultrasoft and hard systems.

We address the frequency dependence of the anharmonicity in the supercooled regime by determining the intrinsic anharmonicity $Q 0$ and compare the experimental findings to MCT predictions of Seyboldt et al. [58]. We extended the experimental results in the supercooled regime to more effective volume fractions and determined the frequency dependence with a higher point density and less scattering compared to the data of Seyboldt et al. [58] by using the frequency sweep MAOS methodology introduced by Singh et al. [70] and could therefore show the overlap of the intrinsic anharmonicity in the high-frequency regime predicted by MCT. Additionally, the MCT predictions [58] of the scaling laws of $Q 0$ in the low-, medium-, and high-frequency regimes agree with our experimental findings, thus increasing the confidence in the linear and nonlinear predictions of the MCT. In addition, we extended the experimental results to the glassy regime.

With this, we provide an extensive dataset on the anharmonicity of a colloidal suspension in both the supercooled and glassy regimes. This dataset will serve as a foundation for future studies using MCT, which are already in preparation.

The most prevalent method for quantifying the frequency-dependent asymptotic nonlinear mechanical behavior of a material is FT rheology, extensively studied by Wilhelm et al. [30,34,39,41]. Note that there are several other methods for analyzing and interpreting the nonlinear oscillatory shear response of a material, as described in Sec. I. This work focuses on FT rheology due to the comparability of the obtained data to the MCT predictions.

The measurement and analysis procedure of FT rheology is shown schematically in Fig. 1(a). In (1) the time-dependent stress signal is shown. By a FT of this stress signal, a magnitude frequency spectrum (2) is obtained. The ratio of the magnitude of the third harmonic to the fundamental $I 3 / 1= I 3/ I 1$ can be used as a measure for the anharmonicity of the material. In the SAOS region, where the anharmonicities are too low to be detected with the specific hardware, a proportionality of $I 3 / 1∝ γ 0 − 1$ is obtained as $I 1∝ γ 0$⁠, and $I 3$ just records the noise level. For higher strain amplitudes, i.e., in the MAOS region also known as the asymptotic region, a proportionality of $I 3 / 1∝ γ 0 2$ is generally expected. For even higher strain amplitudes, i.e., the LAOS region, $I 3 / 1$ levels off to a plateau value. $I 3 / 1$ can be used as a measure of the anharmonicity of the stress response. To quantify the intrinsic anharmonicity of complex fluids as a function of the frequency, a pure frequency-dependent parameter $Q 0(ω)$ introduced by Hyun and Wilhelm [63] can be obtained from the asymptotic regime at low strain amplitudes.

The anharmonic parameter $Q( ω 1 , γ 0)= I 3 / 1 / γ 0 2$ is calculated [see Fig. 1(a4)], enabling the determination of the intrinsic anharmonicity $Q 0$⁠, which is defined as $Q 0( ω 1)= lim γ 0 → 0⁡Q( ω 1 , γ 0)= lim γ 0 → 0⁡( I 3 / 1 / γ 0 2)$ [59]. This procedure is repeated for different frequencies to obtain the frequency dependence of the intrinsic anharmonicity $Q 0( ω 1)$⁠, illustrated in Fig. 1(a5). To obtain the frequency dependence, a strain amplitude sweep of all the relevant angular frequencies $ω 1$ is necessary. Note that we will use the term “strain sweep” instead of “strain amplitude sweep” for brevity in the rest of the work.

Depending on the applied frequency, the investigated amplitude region, and the data point density in the strain sweep, one of these experimental strain sweeps can take up to several days. Therefore, this method is quite time and sample consuming, particularly when exploring a wide frequency range or focusing on the material's low-frequency behavior. For polymer melts, the time-temperature superposition principle is used to reduce the measurement time to obtain the $Q 0$ value over a broad frequency range [59,63]. For colloidal suspensions, the time-temperature superposition principle is not applicable [71]. For water-based systems, the measurement time per sample loading is also limited due to evaporation. Often, evaporation is so pronounced that it is necessary to seal the sample from the environment. If the sample cannot be separated from the sealing material, the sample can only be used within one rheometer loading. To overcome these limitations, the frequency sweep MAOS method was introduced by Singh et al. [70].

In the asymptotic MAOS regime, the anharmonic parameter $Q( ω 1 , γ 0)$ is a plateau value, as shown in Fig. 1(a4). Therefore, rather than extrapolating the amplitude dependence of $I 3 / 1$⁠, the intrinsic anharmonicity can be calculated from one value of $I 3 / 1$ by the following equation:

Consequently, one frequency sweep is sufficient to obtain the frequency dependence of the intrinsic anharmonicity $Q 0$ if the boundaries of the asymptotic region of the sample are known. The boundaries of the MAOS regime can be tested with two strain sweeps at the highest and lowest angular frequency $ω 1$ of interest, as schematically shown in Fig. 1(b).

For Deborah numbers $De= ω 1 τ α$ equal to or greater than 1, the departure from linearity is a $γ 0$-controlled process, whereas for $De≪1$⁠, the departure from linearity is a $γ ˙ 0$-controlled process [72]. Therefore, for $De≈1$ and $De>1$⁠, the strain amplitude of the asymptotic anharmonic behavior is expected to be independent of the frequency $γ 0≠ γ 0( ω 1)$⁠, whereas for lower Deborah numbers $De≪1$ the strain amplitude is expected to be inversely proportional to the frequency $γ 0∝ ω 1 − 1$⁠, as $γ ˙ m a x= γ 0 ω 1$⁠. How our experimental results align with this is discussed in Sec. VI A.

The thermoresponsive polystyrene (PS) poly-N-isopropyl acrylamide (PNIPAm) core-shell particles are synthesized in a two-step synthesis process according to Dingenouts et al. [68]. In the first step, a PS core is synthesized via emulsion polymerization using potassium persulfate (K₂S₂O₈) as an initiator, sodium dodecyl sulfate (SDS) as an emulsifier, styrene as a monomer, and a small amount of N-isopropyl acrylamide (NIPAm) as a comonomer to make the deposition of the oligo-NIPAm in the second synthesis step to the cores more likely. The core is purified by dialyzing against purified water for three weeks. In the second step, a crosslinked PNIPAm shell is synthesized onto the core particles in a seeded emulsion polymerization process. As a crosslinker, $2.5 mol.%$ of N,N′-methylenebisacrylamide (MBA) is used. The initiator K₂S₂O₈ is used as in the core synthesis. The thermoresponsive core-shell particles are purified by ultrafiltration against deionized water for six weeks to remove traces of monomer and free polymer in the suspension. The weighed quantities of the educts are given in Table S1 in the supplementary material.

By this synthesis procedure, two core-shell systems $(V=1.6 L,w≈8 wt.%)$ with different core and shell sizes were synthesized to obtain systems with a certain size ratio. These are used to obtain a bimodal mixture to suppress particle crystallization and to be able to investigate the anharmonic behavior of a glassy system at long time scales. The synthesized systems are listed in Table I.

Particle size distributions were measured by differential centrifugal sedimentation (DCS) on an analytical disc centrifuge (DC24000 UHR, CPS Instruments, Inc., Prairieville, LA, USA). A density gradient was prepared from nine aqueous sucrose solutions (⁠ $1.6 ml$ each) ranging from $8.0$ to $2.0 wt.%$ inside a hollow disk rotating at $24000 rpm$⁠, corresponding to a centrifugal force of 225–277 g acting on the particles. A thin layer of n-dodecane was deposited onto the gradient, thus minimizing evaporation of water, and extending the lifetime of the gradient. The step gradient was allowed to equilibrate within $30 min$⁠, thus yielding a continuous gradient, which is linear in volume. Each particle size measurement was calibrated using polystyrene nanospheres with a nominal size of $251 nm$ and a particle density of $1.054 g c m − 3$⁠. A volume of $100 μ l$ of diluted suspensions of the core-shell particles and their cores $(0.02 wt.%)$ were injected into the spinning disc. During the centrifugation, the particles were sedimenting at rates that depend on the actual size and density of the particles. The measured sedimentation time is inversely related to the sedimentation coefficient [73]. Stokes' law and particle densities measured either on a densitometer or derived from DLS measurements were used to calculate the actual particle sizes from the sedimentation coefficients measured in the DCS experiment. The concentration of the particles arriving at the detector position were determined using light attenuation at $λ=405 nm$⁠, thus allowing for the complete determination of the particle size distribution. Mean particle sizes and dispersities given as weight-average diameter divided by the number-average diameter were obtained this way.

The temperature-dependent DLS measurements were performed with the Zetasizer Nano S from Malvern Panalytical. The wavelength of the laser was $633 nm$ and the scattered light was detected at an angle of $173 °$⁠. Measurements were performed in semi macro cuvettes (PMMA/VWR) in a temperature range from $10$ to $50 ° C$ in steps of $0.5 ° C$⁠. The equilibration time between each step was set to $300 s$⁠. Each measurement was performed three times and the given result is the average of the three measurements. The hydrodynamic radius was determined by cumulant analysis, which is implemented in the instrument software.

The oscillatory and steady shear rheology measurements were conducted on an ARES-G2 rheometer (TA Instruments, Newcastle, USA) using a $40 mm$ plate-plate geometry at a gap of $1 mm$⁠. The samples were put on the plates with a syringe and sealed with silicone oil with a viscosity between 215 and 315 mPa s for the investigated temperatures to suppress water evaporation. All the shear experiments were conducted between $20$ and $10 ° C$ resulting in effective volume fractions between 0.70 and 0.83; see below for its determination. Linear oscillatory shear was measured in an angular frequency range of $ω 1=100 rad s − 1$ to down to $ω 1= 10 − 3 rad s − 1$ depending on the mechanical behavior of the sample at the specific volume fraction and a strain amplitude of $γ 0=1%$⁠. The more solid the samples are, the lower the lowest angular frequency was chosen. All samples were sheared at $γ ˙=100 s − 1$ for $t=2 min$ before all measurements (except for the steady shear measurements) to rejuvenate the sample. The measurements were performed after a waiting time of $t w=90 s$ after preshearing, as within this time the microstructure of the material was rebuilt. For more information, refer to the supplementary material, specifically Fig. S3(c). The steady shear rheology was measured in a shear rate range from $γ ˙=300 s − 1$ to down to $γ ˙= 10 − 4 s − 1$ depending on the rheological behavior of the sample at the certain volume fraction. The LAOS measurements were conducted as strain sweeps from $γ 0=1%$ to $γ 0=5× 10 2%$ for angular frequencies between $ω 1=10 rad s − 1$ to down to $ω 1=3.98× 10 − 3 rad s − 1$ and as frequency sweeps at a fixed strain amplitude of $γ 0=7%$ or $γ 0=10%$ depending on the volume fraction. The frequency sweeps were conducted from $ω 1=10 rad s − 1$ to $ω 1=0.02 rad s − 1$ for all volume fractions. Note that the use of a plate-plate geometry leads to inhomogeneous shear fields, for which corrections are made as described in Sec. VI. Despite this, the plate-plate geometry was chosen because the evaporation of the water could be controlled for longer times as shown in Fig. S3 in the supplementary material.

As a model system, thermoresponsive PS-PNIPAm core-shell particles were synthesized as described in Sec. III A. PNIPAm has a lower critical solution temperature in water, which causes a rapid volume transition at $T L C S T=32 ° C$ of the core-shell systems shown schematically in Fig. 2(a). This reversible volume phase transition is caused by the change in the hydration state of the PNIPAm. By heating, the H-bonds between the water and the amide molecules become weaker, leading to a shrinkage of the crosslinked shell and, at the $T L C S T$⁠, to a complete demixing of the PNIPAm and the water [74], resulting in a collapse of the crosslinked shell around the polystyrene core. This results in a temperature dependence of the hydrodynamic radius $R H$⁠. This temperature dependence measured by DLS is shown in Fig. 2(b) for three different samples: a large particle, a small particle, and a bimodal mixture of the two with a mass ratio of $90:10$ used to suppress particle crystallization. The large and the small core-shell particles have a size ratio of approximately $0.58$ in the investigated temperature range of $T=10−20 ° C$⁠. Between $T=10 ° C$ and $T=25 ° C$⁠, the hydrodynamic radius shows a linear dependence on the temperature with thermal expansion coefficients between $−0.5$ and $−1.0 nm ° C − 1$ depending on the effective crosslinking density followed by a drop with the highest first derivative at $T= T L C S T=32 ° C$⁠. At temperatures above $T>35 ° C$⁠, the hydrodynamic radius $R H$ levels off to a plateau value. From this, it could be concluded that the shell thickness is $29 nm$ for the large and $11 nm$ for the small particles in the shrunken state. To suppress the particle crystallization a bimodal mixture of the large and the small particles is used. The systems used are listed in Table I. All investigations are conducted on the bimodal mixture with a radius ratio of $0.58$⁠. The samples contain potassium chloride $( c K C l=0.05 mol l − 1)$ to screen the electrostatic interactions between the particles. The concentration was chosen according to the literature [75]. It results in a Debye length of $λ D≈1.3 nm$⁠, which is lower than the shell thickness resulting in a purely steric interaction potential. Colloidal stability at this salt concentration is given due to the steric interactions of the particles caused by the polymeric shell.

The linear rheological behavior under oscillatory shear of the bimodal mixture is shown in Fig. 3 for effective volume fractions $ϕ e f f$ below and above the glass transition volume fraction $ϕ g$ for a sample with a weight fraction of $8.3 wt.%$⁠. The volume fraction was varied by changing the radius of the particles by changing the temperature as described in Sec. IV. The measurements were conducted in the linear regime at a strain amplitude of $γ 0=1%$⁠. To allow comparison between the different volume fractions of the sample with differing radii of the colloids, the storage modulus $G ′$ and the loss modulus $G ′ ′$ were normalized to the characteristic energy scale of a particle to the reduced moduli $G r e d ′= G ′( R H 3/ k BT)$ and $G r e d ′ ′= G ′ ′( R H 3/ k BT)$ and the angular frequency $ω 1$ was normalized to the characteristic time scale of the diffusion of a single particle in dilution by its own radius to the rescaled frequency $P e ω=(6π η s R H 3/ k BT) ω 1$⁠. The normalization parameters for the different volume fractions are listed in Table S2 in the supplementary material. It is the same temperature-dependent particle volume $R H 3$ that rescales volume fraction, moduli, frequency, and (see below) shear rate. Note that at high packing fractions, the normalization might not be appropriate, because the hydrodynamic radius was measured with DLS in high dilution. The effective radius of the particles in high concentration may differ due to the finite E-modulus of the hydrogel leading to local deformation of the particles.

In the region between $ϕ e f f=0.75$ and $ϕ e f f=0.78$⁠, measurements were done with small $ϕ e f f$-steps of $Δ ϕ e f f≈0.06$ shown in Fig. S2 in the supplementary material. For the volume fraction of $ϕ=0.76$⁠, no decrease of $G r e d ′$ was found within the investigated time scales. Therefore, the glass transition volume fraction was defined as $ϕ g≈0.76$⁠. This glass transition volume fraction is higher than for monodisperse hard sphere particles, found to be around $ϕ g=0.58$ in computer simulations [77] and experimentally [5,78]. The higher glass transition volume fraction results from the softness and especially the dispersity of the radius of the particles [12,64,65].

The samples with $ϕ e f f> ϕ g$ exhibit a solidlike behavior over the whole investigated frequency range as expected for colloidal glasses [12,64,79], with $G r e d ′= constant$ and $G r e d ′ ′$ being about one decade smaller than $G r e d ′$ and exhibiting a minimum. This minimum is correlated to the transition region between the $α$-relaxation, i.e., the out of cage motion, and the $β$-relaxation of the particles in the glass, i.e., the relaxation or rattling of the particles within the cages [80]. The samples below the glass transition $( ϕ e f f< ϕ g)$ show the same behavior for high rescaled frequencies $P e ω$ with $G r e d ′$ being lower than in the glassy samples. For lower frequencies, (respectively $P e ω$⁠), $G r e d ′$ starts to decrease, while $G r e d ′ ′$ is slightly increasing and then decreasing as well, resulting in a crossover point of $G r e d ′$ and $G r e d ′ ′$ and a maximum in $G r e d ′ ′$ close to the crossover point. Below that, the samples show a flow regime, which is typical for viscoelastic materials. The crossover point is correlated to the inverse of the $α$-relaxation time of the particles in the glass, i.e., the time the particles need to escape their cages. The flow regime shows the possibility of free diffusion of the particles, which escaped their cages.

The $α$- and the $β$-relaxation are getting slower, i.e., shifting to lower frequencies, for increasing volume fraction $ϕ e f f$⁠. For volume fractions higher than the glass transition volume fraction $ϕ e f f> ϕ g$⁠, the $α$-relaxation vanishes completely for our investigation times $(t≈11 h)$ as the particles are densely packed that they cannot escape their cages within the investigation time.

The reduced shear stresses as a function of the shear rate of the bimodal mixture are shown in Fig. 4 for effective volume fractions $ϕ e f f$ below and above the glass transition volume fraction $ϕ g$⁠, determined by oscillatory shear measurements. As for the oscillatory shear measurements, the measured data were double normalized to the radius of the particle, i.e., the stress $σ$ was normalized to the characteristic energy scale of a particle to the reduced stress $σ r e d=σ( R H 3/ k BT)$ and the shear rate $γ ˙$ to the Péclet number $P e 0=(6π η s R H 3/ k BT) γ ˙$ by normalizing to the characteristic time of self-diffusion of a particle in dilute suspension.

For the two lowest investigated volume fractions of $ϕ e f f=0.70$ and $ϕ e f f=0.73$ at low shear rates and consequently low Péclet numbers, the reduced stress is directly proportional to the shear rate, i.e., the viscosity $η$ of the material is constant. The materials exhibit a first Newtonian plateau as the particle's $α$-relaxation time is lower than the inverse of the shear rate leading to free diffusion of the particles. Increasing $P e 0$ leads to a sublinear increase in $σ r e d$ showing that the material is shear thinning. Note that this happens already at $P e 0$ far less than unity. Further increasing the $P e 0$ leads to a steeper increase in $σ r e d$ with approximately $σ r e d∝P e 0 0.5$ predicted for soft particles and associated with interparticle slippage [81–83]. For increasing $ϕ e f f$⁠, the first Newtonian plateau shifts to lower $P e 0$⁠, disappearing from the experimental window for $ϕ e f f>0.75$⁠.

For samples with a volume fraction of $ϕ e f f>0.78$⁠, the reduced stress approaches a constant plateau value at low $P e 0$⁠. At these volume fractions, the $α$-relaxation is slowed down to the point where the particles are unable to diffuse out of their cages and flow without the disturbance due to an external stress of at least the same magnitude as this plateau value. The stress required to cause a material to flow is defined as the yield stress $σ y$⁠. Applying a stress higher than this stress $(σ> σ y)$ leads to a shear thinning of the material.

For higher $P e 0$⁠, the reduced stress $σ r e d$ increases with increasing Péclet number $P e 0$⁠, showing that the shear thinning is getting less pronounced. A yield stress plateau was encountered first at $ϕ e f f=0.78$⁠, close to the experimentally determined glass transition at $ϕ g=0.76$⁠. These flow curves are similar to other investigations on comparable systems of dense colloidal suspensions below and above the glass transition volume fraction [12,79,84].

In Sec. VI, first $I 3 / 1$ in dependence on the strain amplitude for different amplitudes and frequencies is shown for the investigated soft particle system. The behavior of solidlike and liquidlike samples is compared. Additionally, the behavior of the investigated soft system is compared to investigations of Poulos et al. [16] on ultrasoft and hard sphere glass-forming colloids. Then, the frequency dependence of the intrinsic anharmonicity $Q 0$ is investigated for samples below and above the glass transition volume fraction and compared to MCT predictions by Seyboldt et al. [58].

The investigated core-shell system exhibits a liquid to glass transition at an effective volume fraction of $ϕ e f f=0.76$ as shown in Sec. V A. In this section, the LAOS response in the liquidlike and the solidlike states is compared. Note that adjustments have to be made to match the experimental data of plate-plate measurements to cone-plate measurements by accounting for the inhomogeneous flow field in parallel plate geometries [85]. These adjustments can be made by multiplying the strain amplitude of the plate-plate measurements by $0.75$ [34] or multiplying the $I 3 / 1$ values by $1.5$ [85,86] resulting in similar correction factors for $Q 0$⁠. Since this section presents a qualitative comparison between $I 3 / 1$ and the yield strain $γ y$ depending on the state of the sample and the frequency, our direct measurement results are displayed instead of the corrected values. For the quantitative comparison with MCT, adjustments are made according to the description provided later in more detail.

In Figs. 5(a) and 5(b), the linear dynamic frequency sweeps, in Figs. 5(c) and 5(d) the material parameters $G ′$ and $G ′ ′$⁠, and in Figs. 5(e) and 5(f) the relative intensities of the third harmonics $I 3 / 1$ are shown as a function of the strain amplitude $γ 0$ for two different volume fractions in the liquidlike regime $ϕ e f f< ϕ g$⁠. In Fig. 6, the same is shown for two volume fractions in the glassy regime. The two volume fractions below the glass transition show a typical viscoelastic liquid response, while the two volume fractions above the glass transition show the typical characteristics of a soft glass.

For all investigated volume fractions, the relative intensity of the third harmonic $I 3 / 1$ shows a monotonic increase with the strain amplitude where the proportionality $I 3 / 1∝ γ 0 2$ is observed for medium strain amplitudes. This region is defined as the scaling region and is followed by a leveling off to a plateau value.

For the volume fraction of $ϕ e f f=0.71$⁠, most of the investigated frequencies are intentionally chosen to be in the range where the sample shows predominantly viscous behavior. For these angular frequencies, the strain amplitude $γ 0$ of the onset of anharmonicities decreases as the angular frequency $ω 1$ increases, because the maximum shear rate attained during an oscillation cycle is equal to the product of the strain amplitude and the frequency $γ ˙ m a x= γ 0 ω 1$⁠. The anharmonicity levels off to the same plateau value. In contrast, for the two frequencies measured in the predominantly elastic regime of the sample [marked in Fig. 5(a) with the red circle and black square], the plateau value decreases for increasing frequencies. This decrease of $I 3 / 1$ from the plateau value in the LAOS regime is seen for the three volume fractions, measured in the predominantly elastic regime as well [shown in Figs. 5(d), 6(c), and 6(d)] and is in agreement with model predictions [64] and experimental findings on a similar system. It can be rationalized by the transformation of the sample from a viscoelastic solid to a viscoelastic liquid at the yielding point [34].

The plateau value of $I 3 / 1$ is higher for the glassy samples [Figs. 6(e) and 6(f)] than for the liquidlike samples [Figs. 5(e) and 5(f)]. The findings of the soft core-shell systems are in agreement with the findings of Poulos et al. [16] on ultrasoft starlike micelles. It differs from the behavior of hard sphere glasses, where the plateau value and the general course of $I 3 / 1$ show a stronger dependence on the frequency as seen in Fig. 5 of Poulos et al. [16], where the decrease in the plateau value of $I 3 / 1$ is present at low volume fractions for the liquidlike samples and gets even more pronounced in glassy samples [16]. Additionally, the glassy sample in Fig. 5 of Poulos et al. [16] displays a peak in $I 3 / 1$⁠, which has been attributed to the phenomenon of two-step yielding that is observed in hard sphere glasses [15,16].

For the other three volume fractions, the investigated frequencies are chosen to be in the range where the sample exhibits predominantly elastic behavior [see Figs. 5(b), 6(a), and 6(b)]. Here, the onset of anharmonicities (defined by extrapolation of the scaling region of $I 3 / 1∝ γ 0 2$ to $0.005$⁠) shifts to slightly lower strain amplitudes as the frequency decreases, similar to the decrease in the yield strain $γ y$ with decreasing frequency. This will be discussed below. For the two samples in the glassy regime (Fig. 6), the scaling is first $I 3 / 1∝ γ 0 2$ but changes to a scaling with a higher exponent shortly before reaching the plateau for the two lowest frequencies. This transition falls in a similar strain amplitude range as yielding and is also observed in polymeric glasses when the sample deforms plastically and the deformation of the sample becomes localized [87]. In our system, this transition point is shifting to lower strain amplitudes for decreasing frequencies. In Fig. 7(a), the strain amplitude dependence of $I 3 / 1$ for the highest volume fraction is shown at the lowest frequency measured. Here, only a substantial increase in $I 3 / 1$ with the dependence of $I 3 / 1∝ γ 0 4$ is observed as the proportionality of $I 3 / 1∝ γ 0 2$ vanishes at low strain amplitudes into the region where $I 3$ is experimentally too small to be detected above the noise.

In Fig. 7(b), the absolute stress signals of the fundamental and the harmonics are shown for a high volume fraction. The third harmonic of the stress signal shows a power law of $I 3∝ γ 0 5$ instead of the expected power law of $I 3∝ γ 0 3$⁠, while the first and the fifth harmonic still show the expected power laws of $I 1∝ γ 0 1$ and $I 5∝ γ 0 5$ confirming that the unexpected increase in $I 3 / 1$ arises from the increase in $I 3$ rather than from a decrease in $I 1$⁠. This suggests that additional nonlinear processes besides the predicted homogeneous nonlinear viscoelasticity are involved, leading to the observed increase. This amplification might be attributed to the same localization effect of deformation seen in polymeric glasses, where a change in the scaling of $I 3 / 1 ∝ γ 0 2$ to $I 3 / 1∝ γ 0 4$ is observed when the deformation of the sample becomes localized [87].

To the best of our knowledge, an increase with a higher exponent than $I 3 / 1∝ γ 0 2$ was never found before in suspensions. Only recently lower noninteger exponents were found [88,89] and attributed to the particle contact dynamics of attractive suspension networks [89].

The samples in the predominantly elastic regime $( G ′> G ′ ′)$ show an overshoot in $G ′ ′$ as a function of the strain amplitude without showing an overshoot in $G ′$⁠, which is type III behavior according to Hyun et al. [40] and is typical for soft glassy materials as concentrated emulsions [90,91], suspensions [20,92], and colloidal glasses [6,16,64,93]. Recently, Donley et al. [57] could gain more insight into this yielding behavior of type III materials by distinguishing between solidlike and fluidlike contributions to the total energy dissipation combining strain-controlled oscillatory shear with stress-controlled recovery tests. Zero-stress recovery tests [56] allow to distinguish between recoverable and unrecoverable energy dissipation. They could show that the overshoot in $G ′ ′$ comes with a continuous transition from predominantly solidlike to predominantly fluidlike dissipation by showing that the deformations transition from predominantly recoverable to predominantly unrecoverable. In colloidal glasses, a significant recoverable component was observed at flow reversal, at strain amplitudes, where $G ′$ and $G ′ ′$ already crossed [94] in accordance with Donley et al. [95] finding large correlations in the microstructure before and after yielding.

Based on those findings, a microstructural interpretation of the processes during this yielding process can be made. In the low strain region where $G ′$ and $G ′ ′$ are still constant, the strain amplitude is not sufficient to break or irreversibly deform or modify the shape of the local cages of the glassy structure.

The increase in $G ′ ′$ and decrease in $G ′$ for type III materials is directly accompanied by an increase in the unrecoverable strain [57], which characterizes the extent of microscopic changes [95], as the shear amplitude is sufficient to locally deform the cages, but the reformation of the cages is faster than the destruction of the cages. While further increasing the strain amplitude, the destruction of the cages becomes more dominant, leading to a higher unrecoverable strain component than a recoverable component. The maximum in $G ′ ′$ is correlated to the dissipation of energy due to the breaking of the cages for colloidal glasses [64]. This is more pronounced in soft sphere glasses compared to hard sphere glasses due to the deformability and elasticity of the particles.

The peak value of $G ′ ′$ slightly decreases for decreasing frequencies, i.e., the dissipation energy resulting from breaking the cages slightly decreases, as the Brownian motion of the particles supports the opening of the cages [13].

Several methods of the determination of the yield strain $γ y$ from oscillatory shear measurements are frequently used, such as the strain at the crossover point of $G ′$ and $G ′ ′$⁠, at the maximum in $G ′ ′$⁠, at the deviation of $G ′ ′$ from its linear plateau value, or at the intersection between a horizontal line through the plateau modulus with a power-law fit of the high strain behavior of $G ′$⁠. Another method to determine the yield strain is to plot $σ$ vs $γ 0$ on logarithmic coordinates and find the intersection of a line with unit slope at low strains (corresponding to the linear regime) with a power-law fit of the high strain behavior. Those methods were recently compared for different materials [96–99] finding the obtained values for the yield strain and the yield stress to vary about 3 orders of magnitude depending on the method with the yield strain at the crossover of $G ′$ and $G ′ ′$ being the highest. Furthermore, these studies have demonstrated that defining a single yield strain is uncertain. During an oscillation cycle, the material yields and unyields depending on the strain currently applied [97,100].

By looking into the transient behavior during a cycle of oscillation, Donley et al. [97] found the $G t ′$ and $G t ′ ′$ values to show a time dependency as soon as the averaged moduli $G ′$ and $G ′ ′$ show deviations from the linear viscoelastic regime. Additionally, they found that the strain amplitude at which the transient behavior crosses the line of $G t ′= G t ′ ′$ for the first time is close to the strain at the maximum of $G ′ ′$⁠. Therefore, they interpret the strain at which $G ′$ decreases and $G ′ ′$ increases simultaneously as the strain at which small fractions of the material already yield [101], while the material as a whole yields close to the maximum of $G ′ ′$ [97]. Aime et al. [102] found a transition from a slow relaxation mode to a fast relaxation mode in concentrated microgels between the maximum of $G ′ ′$ and the crossover of $G ′$ and $G ′ ′$ associated with yielding. Recently, Kamani et al. [100,103] developed a model with only a single yield stress that can depict the characteristics of the whole amplitude sweep by allowing for a rate-dependent relaxation time and plastic viscosity.

In this work, two different metrics for the yield strain are utilized, namely, the strain at which $G ′ ′$ deviates 5% from its plateau value in the linear regime and the strain at the maximum of $G ′ ′$⁠, respectively. The frequency dependence of the metric $γ y= γ 0( G ′ ′ m a x)$ and the frequency dependence of the metric $γ y= γ 0( G ′ ′ = 1.05 ⋅ G ′ ′ l i n)$ are shown in Figs. 8(a) and 8(b), respectively. Note that both metrics may be based on the same underlying physics as a model with a single yield stress can depict both rheological characteristics [100,103].

For increasing frequency, an increase in the yield strains is observed in accordance with previous experimental findings in soft and hard sphere glasses [10,64,65,104]. This can be rationalized by the cage effect, which is present for parameters, where $G ′> G ′ ′$ in the linear regime. At low applied angular frequencies, the Brownian motion of the particles supports the escape from the cages. At higher frequencies, the particles have less time to escape the cage, i.e., a higher strain amplitude is necessary to break the cage. The suspension behaves still linearly at higher strain amplitudes than it does at lower frequencies.

In Figs. 8(c) and 8(d), the volume fraction dependence of the yield strain amplitudes is shown for three different angular frequencies. It is observed that both yield strains increase with increasing effective volume fraction in accordance with the findings of Koumakis et al. [66] on soft core-shell particles.

In contrast, for hard sphere glasses, a maximum in the yield strain in dependence on the volume fraction was observed determined from oscillatory shear [3,66] and from creep and recovery measurements [65,105]. Note that different methods for the determination of the yield strain $γ y$ do not necessarily yield the same values, but Dinkgreve et al. [96] found the same qualitative behavior for different simple yield stress fluids in dependence on the volume fraction $ϕ e f f$ for different methods. Therefore, a qualitative comparison of the yield strains can be made, but a quantitative comparison between the yield strains obtained with different yielding criteria such as $γ y= γ 0( G ′ ′>1.05⋅ G l i n ′ ′)$ and $γ y= γ 0( G ′ ′ m a x)$ with the values from literature at $γ y= γ 0( G ′ = G ′ ′)$ [3,66] or from creep and recovery measurements [65,105] is not reasonable.

The maximum yield strain observed in hard sphere glasses [65,66] was attributed to two competitive effects. In the glassy regime, the distance between particles decreases while increasing the effective volume fraction approaching random close packing. This results in a decreased strain required to break the cages. Additionally, the yield strain tends toward zero at the liquid-to-glass transition as particles are not permanently caged in the liquid state causing the yield strain to decrease for decreasing volume fraction starting from the maximum.

The increase in the yield strain of soft particle glasses for increasing volume fraction below the glass transition volume fraction can be attributed to the same effect of the emergence of the cage effect as in hard sphere glasses. The increase above the glass transition volume fraction is attributed to the particle softness, which dampens the cage relaxation partially [65,66].

To reduce the measurement time for investigating the frequency dependence of the intrinsic anharmonicity $Q 0$⁠, the measurements can be conducted as a frequency sweep instead of several strain amplitude sweeps as described in Sec. II. To determine the minimum and maximum strain amplitude thresholds for the MAOS region at various frequencies, we utilize the strain sweeps depicted in Figs. 5(e), 5(f), 6(e), and 6(f).

For the effective volume fractions of $ϕ e f f=0.71$ and $ϕ e f f=0.73$ [see Figs. 3, 6(a), and 6(b)], the sample is in the range of $De<1$ for the low frequencies. However, there is still no shift in the strain amplitude range of asymptotic behavior, indicating that the sample is not in the region of $De≪1$⁠. For the higher volume fractions $ϕ e f f=0.80$ and $ϕ e f f=0.83$ in the glassy regime, the sample is in the region of $De>1$ for all investigated angular frequencies.

In the solidlike regime (⁠ $ϕ e f f=0.73$ and $ω 1>0.2 rad s − 1$⁠; $ϕ e f f=0.80$⁠; $ϕ e f f=0.83$⁠), there is a slight shift in the strain amplitude boundaries for asymptotic behavior to lower strains at lower frequencies, as previously described. But it is still evident that for all volume fractions a fixed strain amplitude can be utilized for the frequency sweep MAOS. The strain amplitude was selected between $γ 0=7% and10%$⁠, depending on the volume fraction. Strain sweeps at lower frequencies were conducted within a shortened measurement window to confirm whether the selected strain amplitude is still within the MAOS regime for those low frequencies. Note that for low angular frequencies $( ω 1<0.01 rad s − 1)$⁠, the region where $I 3 / 1$ scales quadratically with $γ 0$ could not be detected as described in Sec. V A with Fig. 7. Therefore, the investigation of $Q 0$ is limited to an angular frequency range of $0.01 rad s − 1≤ ω 1≤10 rad s − 1$⁠.

To further evaluate the suitability of the frequency sweep MAOS and strain amplitude for our system, a comparison between a frequency sweep at a fixed strain amplitude of $γ 0=10%$ and different strain sweeps is shown in Fig. 9 for a volume fraction of $ϕ e f f=0.73$⁠. The $Q 0$ values observed from the frequency sweep and the strain sweeps show the same trend and, especially in the high and low $P e ω$-region, the data overlap. Only in the region around the maximum, the $Q 0$ values determined from the frequency sweep MAOS deviate slightly from the values determined from the strain sweep MAOS. Note that the strain sweep data scatter, while the frequency sweep data are quite smooth. By using a fixed amplitude, the scattering of the $I 3 / 1$ data is reduced, resulting in small deviations from the actual value. However, the values are consistent with each other.

By using this time-efficient method, we could investigate the intrinsic anharmonicity $Q 0$ for a large frequency window and with a high point density in frequency space, in both the liquidlike and the glassy regime while reducing measurement artifacts from water evaporation. Note that determining critical spectra requires wide frequency ranges. In Fig. 10(a), the intrinsic anharmonicity $Q 0$ for the liquid states is plotted against the rescaled frequency $P e ω$⁠. The $Q 0( P e ω)$ curves of the different volume fractions overlap in the high-frequency regime and exhibit a maximum in $Q 0$⁠, which increases and shifts to lower $P e ω$ for higher volume fractions in agreement with MCT predictions from Seyboldt et al. [58] (discussed more in detail below). In addition, predictions from MCT from Seyboldt et al. [58] are shown. The absolute predicted $Q 0$ values based on input parameter from the linear rheology of a disperse system are illustrated as a dashed line, whereas the predicted scaling laws for the low-, medium-, and high-frequency behavior, which are universal for glass-forming colloids, are shown as solid lines.

Due to the inhomogeneity of the flow field in parallel plate geometries, adjustments of the nonlinear rheological results have to be made to match the results to the data from cone-plate and plate-plate experiments. Wilhelm et al. [34] adjusted by multiplying the strain amplitude values of the plate-plate measurements by $0.75$⁠, whereas Wagner et al. [86] and Giacomin et al. [85] multiplied the $I 3 / 1$ values obtained from parallel plate experiments by $1.5$⁠, yielding similar correction factors for Q. We compared the results in $I 3 / 1$ of parallel plate (pp) and cone-plate (cp) experiments and found that in our case, a correction by multiplying the strain amplitude by $0.75$ leads to agreement between the cone-plate and plate-plate results (see Fig. S4 in the supplementary material). A direct multiplication of the $I 3 / 1$ values is not reasonable for our measurements, as the plateau value of $I 3 / 1$ is almost independent of the used geometry. From this, a correction of approx. $Q 0 , C P=1.77⋅ Q 0 , P P$ results. In Fig. 10(b), the adjusted values of $Q 0$ are presented. Additionally, the experimental findings published in Seyboldt et al. [58] are displayed. Agreement between the two datasets is observed for the high-frequency flank with deviations of $6−45%$ between them. Deviations between our work and the work of Seyboldt et al. [58] may be accounted for by the different particle size distributions and modalities, as Seyboldt et al. [58] investigated a polydisperse monomodal system.

The MCT predicts three different regimes with different dependence of $Q 0$ on the rescaled frequency $P e ω$⁠. The predicted scaling $Q 0∝P e ω − a$ with $a=0.32$ in the high-frequency regime and the scaling in the intermediate frequency regime $Q 0∝P e ω b$ with $b=0.69$ agree with our experimental findings. Additionally, a first indication of the predicted quadratic behavior in the low-frequency regime $Q 0∝P e ω 2$ is observed.

Note that the high-frequency regime in $Q 0( P e ω)$ in monodisperse polymer melts shows similar scaling laws and the low-frequency regime shows the same scaling law as in colloidal suspensions, so an exponent of $0.35$ and $2$⁠. In polymer melts, the exponent of the low-frequency regime highly depends on the dispersity of the molecular weight of the system [59].

The different regimes shift to lower angular frequencies, respectively, lower $P e ω$⁠, for higher volume fractions consistent with MCT predictions as the relaxation times of the sample are increasing for increasing volume fractions.

For high $P e ω$⁠, the MCT theory predicts the overlap of the intrinsic anharmonicity $Q 0$ of the different volume fractions. Our findings validate this overlap and coincide with the findings of Seyboldt et al. [58] on a similar system. Note that the data shown in Seyboldt et al. [58] show the expected proportionality in the high-frequency region only for the highest investigated volume fraction $ϕ e f f$⁠, while our experimental data show the dependence for all investigated volume fractions as our experimental results show reduced scattering compared to the data of Seyboldt et al. [58] due to the use of the time-efficient frequency sweep MAOS instead of using strain sweep MAOS. The deviations of the experimental results from the absolute predicted values by MCT [depicted as a dashed line in Figs. 10(a) and 10(b)] can be explained by the difference in the particle size distribution, as Seyboldt et al. [58] investigated a monomodal, but disperse suspension $( σ r e l=17%)$⁠, i.e., the input parameter of the MCT predictions was based on the linear rheology of this system. In contrast, the experimentally investigated system in this work is bimodal.

The similarities between the experimental findings and the MCT predictions imply the need for further scrutiny of MCT predictions, which are based on input parameters from the linear rheology of the bimodal system. In addition, the qualitative agreement between the MCT predictions, which were often used to describe hard sphere glass-forming colloids, and our experimental results of soft core-shell systems further supports the independence of the anharmonicity from the particle interactions, which was proposed by Poulos et al. [16]. They found qualitative agreement of the anharmonicity of the two more extreme cases of ultrasoft starlike micelles and hard sphere glasses.

In Fig. 11(a), the course of $Q 0$ for the liquidlike and the glassy samples is shown as determined from the frequency sweep MAOS via $Q 0= I 3 / 1 / γ 0 2$⁠. The glassy samples show the predicted scaling $Q 0∝P e ω − a$ with the critical exponent $a=0.32$ in the high-frequency regime in quantitative agreement with the liquidlike samples. The intrinsic anharmonicity $Q 0$ of the liquidlike samples overlap, whereas in the glassy regime, $Q 0$ decreases for increasing volume fractions, which is not expected. Therefore, in Fig. 11(b) the intrinsic third anharmonicities $[ I 3]= I 3 / γ 0 1$ and the intrinsic fundamentals $[ I 1]=|G∗|= I 1 / γ 0 1$ are shown to analyze if the decrease arises from a decrease in $[ I 3]$ or an increase in $[ I 1]$⁠. They can be used as an internal comparison, since the samples are measured with the same rheometer, and thus with the same torque sensitivity.

We can conclude that the decrease in $Q 0$ is due to the increase in $[ I 1]$ rather than a decrease in the intrinsic third harmonic $[ I 3]$⁠. $[ I 3]$ and $[ I 1]$ increase with the volume fraction $ϕ e f f$ in the liquidlike region, whereas in the glass this increase vanishes for $[ I 3]$ while $[ I 1]$ still scales with the volume fraction $ϕ e f f$⁠. Consequently, the decrease in the intrinsic anharmonicity $Q 0$ is caused by the higher linear elasticity of the sample approaching higher volume fractions rather than a decrease in the intrinsic third harmonic $[ I 3]$⁠. The collapse of $[ I 3]= I 3/ γ 0 3$ in the glass is suggestive and should be analyzed by MCT predictions based on the linear rheology of our sample.

In Fig. 11(d), the volume fraction dependence of $P e ω$ of the characteristic linear and nonlinear mechanical properties is shown. The minimum in the loss factor $tan⁡(δ)$ and in $G ′ ′$ extracted from the linear frequency sweeps (see Sec. V A, Fig. 3) shows different scaling laws in the liquid and in the glass. The maximum in $Q 0$ obtained from the frequency sweep MAOS measurements [see Fig. 10(a)] lies in the region of the $α$-relaxation process and scales accordingly. The maximum in $Q 0$ is about two decades higher in amplitude than that found in linear and sparsely branched polymer melts [59]. This is in agreement with the findings on soft colloidal star polymer systems [106] and is attributed to the difference in the onset of nonlinearities. Polymer melts behave linearly over a wider range of strain amplitudes $( γ y≈100%)$ than colloidal suspensions $( γ y≈10%)$ resulting in lower maximal $Q 0$ values compared to colloidal suspensions. The $Q 0$ value is lower for a system that exhibits linear behavior over a wide range of strain amplitudes compared to a system, where the onset of nonlinearities occurs at smaller strain amplitudes.

In agreement with Poulos et al. [16], we find an increase in $I 3 / 1$ for samples in the supercooled state in the frequency region, where the sample shows a predominantly viscous behavior in the linear regime and a decrease for frequencies in the region, where the sample shows predominantly elastic behavior in the linear regime. This can be explained by the increasing shear rates reached during a cycle of oscillation. For liquidlike samples, this leads to stronger shear thinning, resulting in an increase in the nonlinearities. In contrast, for solidlike samples, this leads to a transition from elastic to plastic flow and, therefore, to a decrease in anharmonicities as the material cannot follow the applied external shear rate [16]. The intrinsic anharmonicity $Q 0$ is closely related to the linear frequency-dependent behavior, as it describes the anharmonicity of the sample in the limit of vanishing strain amplitude $Q 0= lim γ 0 → 0⁡( I 3 / 1/ γ 0 2)$⁠. For volume fractions in the glassy state, no maximum in $Q 0$ is observed supporting the MCT predictions of a power-law divergence for low frequencies at the glass transition. This demonstrates the need for MCT predictions of the glassy state.

Soft thermoresponsive PS-PNIPAm core-shell particles in suspension were synthesized and purified to investigate the nonlinear oscillatory rheology of the glass-forming suspension below and above the glass transition volume fraction for answering a question, which was raised by previous comparison of FT rheological results to MCT [58]. It addresses the behavior in supercooled states close to the glass states for which we provide data in a broad frequency range. To suppress particle crystallization, a bimodal mixture with a size ratio of approx. $0.58$ in radii in the investigated temperature range of $10−20 ° C$ and a mass ratio of large $( R H , 20 ° C=135 nm)$ to small $( R H , 20 ° C=79 nm)$ particles of $90:10$ was used, instead of using a disperse system enabling better comparability to MCT predictions. The nonlinear rheology of the samples was investigated at both low and high volume fractions, where the behavior is liquidlike and glassy, respectively, using frequency and strain sweep MAOS.

The soft core-shell particles behave similarly to the ultrasoft system studied by Poulos et al. [16]. Their study investigated the nonlinear behavior of ultrasoft starlike micelles and hard sphere glasses. We found that the core-shell particles feature a constant plateau value of $I 3 / 1$ for angular frequencies in the predominantly viscous regime. For angular frequencies in the dominantly elastic regime, there is a slight decrease in this plateau value with increasing $γ 0$⁠. The ultrasoft systems show the same behavior. In contrast, the hard sphere glass-forming colloids already show a decrease in the plateau value in the predominantly viscous regime and a significant drop in the glassy regime [16].

In addition, a substantial increase of $I 3 / 1$ at low angular frequencies with a dependence of $I 3 / 1∝ γ 0 4$ in colloidal glasses was observed. This behavior is similar to that observed in amorphous polymeric materials when they undergo plastic deformation, where the deformation of the sample becomes localized [87].

The oscillatory yield strain amplitude (obtained from the strain at the maximum in $G ′ ′$⁠) increases for increasing frequencies and for increasing volume fraction. The increase in yield strain for increasing frequencies can be attributed to the reduced time the particles have to escape from their cages. Thus, higher strain amplitudes are required to break the cages.

Additionally, we supported the findings of Singh et al. [70] that utilizing frequency sweep MAOS is a viable and more time-efficient alternative to strain sweep MAOS to obtain the frequency-dependent anharmonic characteristics of colloidal glasses. By using frequency sweep MAOS instead of strain sweep MAOS, we could determine the anharmonicity of glass-forming colloidal suspensions in a broad frequency and volume fraction regime with a high point density in frequency space, which makes a thorough comparison to MCT predictions possible. Our experimental results show little deviation (max. $∼50%$⁠) from the quantitative MCT predictions of $Q 0$ in the high-frequency region. This supports the MCT prediction of a power-law rise of the nonlinearity with decreasing frequency that is independent of the distance to the glass transition.

Furthermore, our findings reveal that $Q 0$ exhibits the scaling laws predicted by MCT in the three different frequency regions [58]. The frequencies at which these scaling laws occur depend on the radius and volume fraction of the particles. The qualitative agreement between the MCT predictions of $Q 0$ based on a monomodal, disperse system with the experimental results of the bimodal system investigated in this work shows the universality of the critical nonlinearity at the MCT glass transition and suggests a more detailed analysis of MCT predictions, which start from the measured linear rheology of the bimodal system.

It was demonstrated that for colloidal liquids $Q 0$ reaches a maximum in the region of the $α$-relaxation, and its position scales accordingly. As the sample becomes more fluidized, the material can dissipate the energy introduced via the shear, which leads to a decrease in the intrinsic anharmonicity $Q 0$⁠. In contrast, in the solidlike regime, the higher shear rates reached during a cycle of oscillation lead to a transition from predominantly elastic to predominantly viscous flow, which leads to a monotonous decrease in the anharmonicity.

The glassy samples still show the predicted scaling of $Q 0∝P e ω − 0.32$⁠. In the high-frequency limit, the intrinsic anharmonicity $Q 0$ is not a function of the volume fraction in the liquidlike region, whereas it decreases with volume fraction in the glassy region. In contrast, the intrinsic third harmonic of the stress signal $[ I 3]$ increases in the liquidlike region, whereas above the glass transition, it is no longer a function of the volume fraction, while the intrinsic fundamental $[ I 1]$ and consequently the stress $σ$ increases with increasing volume fraction. Consequently, the decrease in the intrinsic anharmonicity $Q 0$ is caused by the higher elasticity of the sample approaching higher volume fractions rather than a decrease in the absolute third harmonic $I 3$⁠. The collapse of $[ I 3]$ in the glass indicates the need for a comprehensive analysis of MCT predictions based on the linear rheology of our sample.

See supplementary material for the synthesis details; normalization factors for the frequency, shear rate, moduli, and stress; frequency sweeps with small effective volume fraction steps to investigate the glass transition volume fraction with high precision; time sweeps of the sample with plate-plate and cone-plate geometries to show the sample stability over time and the time scale of the recovery of the microstructure and the comparison between the relative intensity of the third harmonic measured with a plate-plate and with a cone-plate geometry.

The authors gratefully acknowledge the financial support of the German Research Foundation (DFG) through Project Nos. WI 1911/28-1 and FU 309/12-1. Furthermore, the authors thank Matthias Karg and Marco Hildebrandt for conducting the temperature-dependent DLS measurements and Michael Pollard for proofreading the manuscript as a native English speaker.

The authors have no conflicts to disclose.

The data that support the findings of this study are available from the corresponding author upon reasonable request.