Predicting nonlinear physical aging of glasses from equilibrium relaxation via the material time

The noncrystalline glassy state of matter plays a role in virtually all fields of materials science and offers complementary properties to those of the crystalline counterpart. The caveat of the glassy state is that it is out of equilibrium and therefore exhibits physical aging, i.e., material properties change over time. For half a century, the physical aging of glasses has been known to be described well by the material-time concept, although the existence of a material time has never been directly validated. We do this here by successfully predicting the aging of the molecular glass 4-vinyl-1,3-dioxolan-2-one from its linear relaxation behavior. This establishes the defining property of the material time. Via the fluctuation-dissipation theorem, our results imply that physical aging can be predicted from thermal-equilibrium fluctuation data, which is confirmed by computer simulations of a binary liquid mixture.


INTRODUCTION
Physical aging deals with small property changes resulting from molecular rearrangements [1][2][3][4][5].While the aging of a material is, in practice, often due to chemical degradation, physical aging does not involve any chemical change.Understanding this type of aging is crucial for applications of noncrystalline materials such as oxide glasses [4,[6][7][8], polymers [2,5,[9][10][11][12][13], metallic glasses [14][15][16][17][18], amorphous pharmaceuticals [19], colloidal suspensions [20], etc.For instance, the performance of a smartphone display glass substrate is controlled by details of the physical aging during production [21], and some plastics eventually become brittle as a result of physical aging [22].Noncrystalline or partly noncrystalline states play a role in modern materials science, e.g., in connection with metal-organic frameworks [23] and high-entropy alloys [24], and physical aging is also important in connection with active matter [25][26][27].Last, it should be mentioned that aging under nanoconfined conditions differs from that of bulk materials [28].The lack of a fundamental understanding of the glassy state and its aging with time influences all branches of materials science, which explains the continued interest in the field from a theoretical point of view [17,27,[29][30][31][32].
Describing and predicting physical aging has been a focus of glass science for many years, yet the subject still presents important challenges [13,31].In this work, we address the concept of a material ("reduced") time controlling aging, which was proposed by Narayanaswamy in 1971 in a paper dealing with the physical aging of oxide glasses [6].A closely related formalism describing polymer aging was developed a few years later by Kovacs and coworkers [3] and, in the 1990s, in the entirely different context of spin glasses by Cugliandolo and Kurchan [29].The materialtime concept rationalizes several notable aging phenomena [4,6,[33][34][35].It is used routinely in both basic research and applications.The material-time formalism is generally recognized to describe well the physical aging of systems subjected to relatively small temperature variations, but the existence of a material time has never been validated in direct experiments.We do this here in long-time experiments on a glass-forming molecular liquid by demonstrating the fundamental prediction that linear response aging data determine the nonlinear aging behavior in the intermediate regime involving temperature variations of a few percent.
Physical aging is a complex phenomenon as it is both nonexponential and nonlinear.The simplest and best controlled aging experiment is based on the temperature jump protocol: The sample is initially in a state of thermal equilibrium, then its temperature is changed instantaneously, i.e., rapidly compared to the response time scale of the material, and the full approach to equilibrium at the new temperature is monitored as a function of time [36].This procedure requires a setup that allows for fast temperature changes and has a precise temperature control with a minimal long-time drift.Moreover, accurate measurements are needed because the long time tail of physical aging, as well as the entire aging response to a small temperature step, involves only minute changes of material properties.
Our experimental setup is based on a Peltier element in direct contact with a plane-plate capacitor.The setup keeps temperature constant over months with less than 1 mK variation, and the samples are so thin (50 µm) that the temperature may be changed within a few seconds to a new, constant value.Dielectric properties are monitored using an ultra-precision Andeen-Hagerling capacitance bridge.More details on the setup are provided in Materials and Methods and in 36-40.

RESULTS
We performed several temperature jump experiments around a reference temperature on the glass-forming liquid 4-vinyl-1,3-dioxolan-2-one (VEC) and monitored, after each jump, both the real and the imaginary part of the capacitance at 10 kHz as the system gradually equilibrates [41,42].The real part of the VEC data is presented in Fig. 1; the imaginary part of the data can be found in the Supplementary Materials in which we also give analogous data for N-methyl--caprolactam (NMEC).Capacitance can be measured very precisely and is an excellent probe in aging experiments [36,43,44].For samples of molecules with a low dipole moment, the real part of the capacitance provides a direct measure of the density [38].The VEC and NMEC molecules have large dipole moments, which implies that rotational polarizations contribute substantially to the capacitance even at high frequencies [38,45].
Although this means that the simple connection to density is lost, the capacitance still provides a precise probe of the state of the sample during aging.The reference temperature for the VEC experiment is 164.6 K at which the main (alpha) relaxation time is roughly 12 hours (see the Supplementary Materials).This is large enough for the setup to thermalize after a temperature jump before any significant relaxation has taken place in the sample.Figure 1A shows our temperature protocol with the 10 kHz real part of the capacitance measured as a function of time.The first 36 weeks of the experiment were devoted to single temperature jumps with size varying from 10 mK to 3 K, carefully equilibrating the sample after each jump before the next one was initiated.The last 15 weeks were spent on temperature variations involving double jumps and sinusoidal modulations.The latter are not resolved in this figure, where they appear as thick turquoise vertical lines; we return to these protocols later (Figs. 2 and 3). Figure 1B shows the data for the single jumps plotted as a function of the logarithm of the time that has passed after each jump was initiated.Note that these curves have very different shapes, demonstrating that even fairly small temperature jumps lead to a notable nonlinear response.This is a hallmark of physical aging, reflecting the "asymmetry of approach" that jumping to the same final temperature from a higher temperature results in a faster and more stretched response than the same size jump coming from below [4,6,13,46].Figure 1C focuses on the smaller jumps that are not resolved in Fig. 1B.
The response to a temperature variation is usually highly nonlinear.Nevertheless, any response is expected to have a small-amplitude limit at which the measured quantity, X(t), depends linearly on the external perturbation.Even in this limit, the measured signal, in general, depends on the temperature history.This means that X(t) in the linear limit is given by a convolution of the change in temperature with the normalized linear time-domain response function R lin (t) in the following manner Here, the temperature T is in general a function of the time t, X eq (T ) is the equilibrium value of the measured property at temperature T , and α X = dX eq /dT quantifies its temperature dependence.This linear description is also known as the Boltzmann superposition principle.In the case of a temperature jump at time zero from the initial temperature T i to the "bath" temperature T b , one has dT /dt = −∆T δ(t ) in which (following the convention in the field) ∆T = T i − T b and δ(t ) is the Dirac delta function.Equation 1 implies that the time-dependent response is α X ∆T R lin (t).Defining ∆X = X eq (T i ) − X eq (T b ) and noting that ∆X = α X ∆T , the response is given by We have previously worked with this linear limit for temperature jumps down to 100 mK [38][39][40]47]; the data of the present paper take this a step further by involving temperature jumps as small as 10 mK, as well as by optimizing the protocol to make it possible to properly resolve both the long-and the short-time plateaus of the linear aging curve.Linearity is investigated, in general, by considering the normalized relaxation function of the quantity X, denoted by R(t), which for a jump to temperature T b at t = 0 is defined by R(0) = 1 and R(t) goes to zero as the system equilibrates at T b at long times.Whenever the data are in the linear regime, the relaxation function is the response function of Eq. 1, R(t) = R lin (t), i.e., relaxations following all temperature jumps have the same time-dependent normalized relaxation function in the linear limit.
Figure 2 shows R(t) for temperature jumps of magnitude 10 to 100 mK to and from the reference temperature, where X = C (10 kHz) is the real part of the capacitance at 10 kHz of VEC.Similar data are shown for the imaginary part and for NMEC in the Supplementary Materials.The short-time plateau of R(t) is below the theoretical value R(0) = 1.This is because there is a fast response that cannot be resolved by our setup, a common finding in studies of physical aging that was discussed in detail in previous works [39,47].We believe that the fast response mainly happens on the phonon time scale due to vibrational and librational equilibration.In addition, one or more beta relaxations may take place at times shorter than the experimental cutoff of about 4s (the time it takes to change temperature and equilibrate the setup at the new temperature).However, in the dielectric spectra (fig.S1), the beta relaxation is only seen as a small shoulder, which indicates that it probably only accounts for a few percent of the initial decay in R(t).All data collapse as predicted for linear relaxation.B Data from a small-amplitude temperature double jump starting at T0 = 164.6K and jumping first by -100 mK and then by +50 mK (colored curves) along with the prediction according to Eq. 4 (black dashed lines).The upper panel shows the full experiment on a linear time scale, the lower panel shows the data on a logarithmic time scale that sets the time of the beginning of each temperature jump to zero.The inset illustrates the temperature protocol.C Sinusoidal small-amplitude temperature protocol (red points) and data (turquoise points behind the black points).The prediction based on the linear-response formalism (Eq.5) is shown as black points.The deviations between prediction and data that can barely be discerned in the main figure can be seen in the magnification.
All the normalized relaxation functions observed for these small temperature jumps around the reference temperature collapse within the experimental uncertainty, as predicted for a linear response.True linearity is a theoretical limit, which means that higher precision and better resolution would reveal tiny differences between the relaxation curves.For the data of Fig. 2, the uncertainty is of the same order of magnitude as the symbol size and no differences are resolved, meaning that the measured curves for all practical purposes represent the linear-response function R lin (t).In the following, we demonstrate how the linear response function can be used to predict the response for different temperature protocols resulting in linear (Fig. 2) and nonlinear (Fig. 3) aging responses.
A simple generalization of the temperature-jump experiment is to introduce a second jump before the system has equilibrated fully in response to the first one, a so-called double-jump experiment.If the temperature changes are both small enough to be within the linear range, then the measured value of X(t) after the second jump is a sum of the responses to the individual jumps.For two temperature jumps corresponding to changes in the measured property X by ∆X 1 and ∆X 2 performed at times t 1 and t 2 (t 1 < t 2 ), respectively, one has where T 2 is final temperature after the second jump.We test Eq. 4 for the Ritland-Kovacs crossover protocol [4,46,48,49] consisting of two consecutive temperature jumps with an opposite sign determined such that the observable X has its equilibrium value right after the second jump.Figure 2B illustrates this protocol and shows the observations after a −100 mK jump followed by a +50 mK jump for VEC.The data reproduce the crossover effect that X(t) exhibits a peak after the second jump [4,46].This bump is a manifestation of the memory present for any nonexponential linear response [4].Along with the data, the predictions based on Eq. 4 and the measured R lin (t) from the 50 mK jump in Fig. 2A are also shown.The prediction collapses almost exactly with the double-jump data.These small-amplitude double jump results provide an extra confirmation that we have reached the linear limit of physical aging.Similar data are presented for NMEC in the Supplementary Materials, which also provides data for more VEC small (linear) jumps.
Moving on from the double temperature jump, we note that Eq. 1 predicts the response to any temperature perturbation small enough to be linear.Because we do not have an analytic expression for R lin (t), the integral is calculated by generalizing the sum in Eq. 4: where T N is the final temperature after N jumps.In Fig. 2C, we show how Eq. 5 predicts the output of a small sinusoidal temperature perturbation.The frequency of the perturbation is 2.3 • 10 −5 Hz, which is the inverse of the estimated equilibrium relaxation time of the sample at the reference temperature 164.6 K.The amplitude is 100 mK, i.e., within the linear regime of single jumps.The prediction follows the data with a high accuracy, including both the transient behavior (seen, e.g., in a first peak that is higher than the second) and the phase shift.Tiny deviations between prediction and data can be seen in the inset, which also shows how the temperature protocol is composed of 2 mK temperature steps.After establishing the linear aging limit and showing how linear temperature-jump data can be used to predict the response of other linear temperature protocols, we now turn to the main result of this paper, a proof of the existence of a material time for VEC.The radically new idea in the 1970s [3,6] was that aging becomes linear when it is described in terms of the material time ξ(t) instead of the laboratory time t.One assumes the so-called time aging-time superposition, meaning that the spectral shape of R lin is independent of the state of the sample.As a consequence of these assumptions, Eq. 1, Eq. 4, and Eq. 5 describe also nonlinear experiments by replacing the laboratory time with the material time, i.e., The material time is "measured" by a clock with a rate that reflects the state of the sample, and the nonlinearity of physical aging is a consequence of this fact [3,6,29,50].The material time may be thought of as analogous to the proper time in the theory of relativity, which is the time recorded on a clock following the observer.Although a microscopic definition of the material time remains elusive, this concept is generally recognized to form the basis of a good description of physical aging involving relatively small temperature variations [4].The very fundamental assumption of the formalism, however, that nonlinear aging phenomena can be predicted from the linear aging limit, has never been validated.In the following, we do so by showing how the measured linear response determines the response to nonlinear temperature protocols for VEC and, in the Supplementary Materials, for NMEC.
3. Data from the real part of the capacitance, C (10 kHz), from large-amplitude temperature-modulation experiments on VEC along with predictions based on the measured linear response R lin (t), i.e., using Eq.6 in combination with Eq. 7 and Eq. 8. Colored curves are data and black dashed lines are predictions.A Normalized relaxation functions (Eq. 3) of single temperature jumps with amplitudes ranging from 1 K to 3 K.B Data from a double jump starting at T0 = 167.6K and jumping by -2 K and +1 K (colored curves).The upper panel shows the full experiment on a linear time scale, the lower panel shows the data on a logarithmic time scale, setting the time of the beginning of each temperature jump to zero.The inset illustrates the temperature protocol.C Sinusoidal temperature protocol (red points) and data (turquoise points) shown together.The horizontal dashed line marks the equilibrium capacitance at the starting temperature 164.6 K.The prediction of the material-time formalism is shown as black points.The response is highly nonlinear, resulting in a non-sinusoidal curve that is far from symmetric around the horizontal dashed line.The deviations between prediction and data that can barely be discerned in the main panel are visible in the magnification.
Using Eq. 6 requires a connection between the laboratory time t and the material time ξ.This is obtained by introducing the time-dependent aging rate γ(t) defined [3-6, 8, 11, 51] by In equilibrium, the aging rate equals the relaxation rate γ eq defined as the inverse of the equilibrium relaxation time.Thus, a linear experiment is the limiting case for which the aging rate is constant and the material time is proportional to the laboratory time, ξ lin (t) = γ eq t.
Here, γ eq (T ) and X eq (T ) are the equilibrium values of γ and X at the temperature T , and Λ is a constant that depends only on the substance and the monitored property X.It should be noted that Eq. 8 is arrived at by first-order Taylor expansions and, for this reason, can only be expected to apply for relatively small temperature variations.
The material-time description in Eq. 6, combined with Eq. 7 and Eq. 8, gives a unique prediction for X(t) for any temperature protocol.Equation 6predicts the value of X(ξ) while Eqs.7 and 8 connect the material and laboratory times by stretching or compressing the time-scale axis.The input needed for the prediction is R lin (t) as determined in Fig. 2A, the equilibrium values of the rate γ eq (T ) and of the measured property X eq (T ), and the parameter Λ.We have equilibrium measurements of X eq (T ) down to 163.6 K and have extrapolated values to lower temperatures (see the Supplementary Materials).The values used for γ eq (T ) are extrapolations from a fit of relaxation times derived from dielectric spectra, which, down to 163.6 K, are proportional to the aging rates (see the Supplementary Materials).
The parameter Λ is determined by the method described in Ref. 53 from the two temperature-jump experiments of magnitude ±1 K to the reference temperature 164.6 K (see the Supplementary Materials).This Λ value was used for predicting all other nonlinear responses.Figure 3A shows R(t) data from the nonlinear single temperature jumps.It is seen that the short-time plateaus of R(t) for the different jumps do not coincide.This is due to a difference in the short-time relaxation deriving from the response on the phonon time scale and, possibly, also from one or more beta relaxations.To predict the aging, we have adjusted for this difference in a manner where the short-time decay of R(t) depends on both the initial and final temperatures (see the Supplementary Materials).
Figure 3 reports the main results of the paper: data from nonlinear temperature protocols along with predictions based on the linear temperature-jump data.The nonlinear protocols mirror the linear protocols of Fig. 2. Figure 3A shows single temperature jumps, Fig. 3B shows a 2 K and 1 K double jump, and Fig. 3C shows a sinusoidal temperature modulation with amplitude 3 K and the same frequency as the linear sinusoidal protocol of Fig. 2.
The single jumps in Fig. 3A exhibit the asymmetry of approach characteristic of nonlinear aging [13,46]: "selfacceleration" of up jumps where the relaxation rate speeds up as equilibrium is approached and "self-retardation" of down jumps [11].The material-time formalism captures well this asymmetry (black dashed lines), and the measured data are predicted with a high accuracy for all down jumps and for up jumps up to 2 K.However, there is a clearly visible deviation for the largest (3 K) up jump and in the Supplementary Materials it is documented that deviations in fact emerge already for a 2.5 K up jump.Thus, the formalism breaks down for large amplitude up jumps.This may be related to only going to first order in the Taylor expansion in Eq. 8, but it could also be caused by the sample reaching equilibrium by other mechanisms than the one involved in smaller jumps.This may be similar to what is seen in the case of very large up jumps (30 to 70 K) performed on ultrastable vapor-deposited glasses where it has been shown that equilibrium is reached by heterogeneous growth of mobile domains [55].Alternatively, the deviations between data and predictions could be caused by beta processes playing a role in aging, as has been seen for polymers deep in the glass state [56,57].
The predictions agree well with the data of the nonlinear double jumps shown in Fig. 3B.This demonstrates that the material-time formalism works well also in this situation; we note that Ref. 41 presents an alternative approach for predicting the nonlinear aging response from linear data.The data shown in Fig. 3B are all from measurements in the temperature range where we have access to measured values of X eq and to the fast contribution of R(t), while γ eq used for the prediction is derived from an extrapolation of higher-temperature dielectric relaxation times.The parameter Λ is the same as for the single jumps, and the test of the nonlinear double-jump prediction is therefore performed with no free parameters.In contrast, in the classical Ritland-Kovacs crossover experiment [46,48,49], the first down jump goes deep into the glass state where the properties of the equilibrium liquid are not known.In the Supplementary Materials we show data for a large down jump (7 K); the predictions using extrapolated parameters demonstrate qualitatively good results, although the formalism is not able to predict the time scale of aging in the temperature regime where equilibrium cannot be reached.
Last, Fig. 3C shows the response of the nonlinear sinusoidal temperature modulation along with the predictions.The lowest temperatures in the modulation are in a range where the parameters are extrapolated.Again there are no free parameters in the prediction.The nonlinearity is seen as a sizable asymmetry in the peak shape: When the temperature is high, there is a substntial response whereas the liquid responds much less to a decreased temperature.The gray horizontal dashed line corresponds to the equilibrium capacitance at the starting temperature 164.6 K.The asymmetry of the response is very well captured by the prediction.Because a large part of a sinusoidal is close to linear in time, and thus similar to a temperature ramp over several Kelvin, the aging in connection with a standard differential scanning calorimetry cooling or heating protocol is likewise expected to be predicted accurately.Deviations between prediction and data can be seen in the magnifications of Fig. 3C and are most likely related to the first-order nature of Eq. 8.The results in Fig. 3 demonstrate that nonlinear physical-aging phenomena in the intermediate regime may be predicted from a knowledge of the linear limit of aging.Previous works have come close to this limit [39,58].While the linear limit is challenging to probe experimentally, it is conceptually important.First, it validates the central assumption of the material-time formalism.Second, the linear response theory is well established via the fluctuation-dissipation (FD) theorem that predicts the response from thermal equilibrium fluctuations quantified via a time-correlation function [59].Our results therefore imply that intermediate nonlinear physical aging can now, at least in principle and for relatively small jumps, be predicted from measurements of the equilibrium fluctuations, i.e., without perturbing the system at all.We end the paper by illustrating this possibility by presenting results from a computer simulation where thermal fluctuations are much easier to monitor than in experiments.
The system studied is the binary Lennard-Jones(LJ) mixture of Kob and Andersen [60], which, for more than 20 years, has been the standard model for computer simulations of glass-forming liquids.We simulated a system of 8000 particles.The quantity monitored is the potential energy U .Temperature-jump data were averaged over 1000 simulations to reduce the noise.Figure 4A shows results for jumps from four different temperatures to T = 0.60 (in simulation units), plotted as a function of the logarithm of the time passed after each jump was initiated.The curves are quite different, showing that the jumps are large enough to be notably nonlinear.
The FD theorem implies that the linear response to any small temperature variation is uniquely determined by the thermal-equilibrium potential-energy time-autocorrelation function U (0)U (t) [61].We evaluated this quantity at T = 0.60.Using the single-parameter material-time formalism as above we then predict nonlinear temperature-jump results (Fig. 4B).The only free parameter is the Λ of Eq. 8, which is determined from the two smallest jumps [53].The colored dashed curves in Fig. 4B are the predictions for the normalized relaxation functions based on the black line in the middle that gives the thermal equilibrium normalized time autocorrelation function of the potential energy; the full circles are the normalized data from Fig. 4A.Overall, the predictions work well, demonstrating that intermediate nonlinear aging can be predicted from equilibrium fluctuations.The minor deviations for the two largest jumps are not unexpected, given that these involve temperature changes of more than 15% for which the single-parameter ansatz in Eq. 8 is likely not to be accurate.

DISCUSSION
We have shown how physical aging involving temperature changes of a few percent can be predicted from the linear aging response, i.e., from the response to a very small temperature variation.This validates the central assumption of the material-time formalism.At the same time, it is clear that this formalism has limitations.Thus, the largest up jump (3 K) is not well predicted (Fig. 3A, bottom).This suggests that there are two regimes of nonlinear aging: an intermediate regime where the relaxation time variesfor, at most, a few decades and the material-time concept describes the situation well, and a strongly nonlinear regime where the formalism breaks down and a new theoretical approach is needed.We speculate that even very large temperature down jumps may fall into the intermediate regime because the system here thermalizes gradually.This is in contrast to large up jumps, which are known to result in heterogeneous states very far from equilibrium [55].Aging far below the glass transition is also likely to deviate from the predictions because processes faster than the alpha relaxation may play a role here, particularly for polymers [56,57].Along this line of thinking, it is important to note that the standard glass transition resulting from a continuous cooling is likely to be described well by the material-time formalism, i.e., is intermediately nonlinear because vitrification for a constant cooling rate takes place over a narrow range of temperatures.
In regard to the intermediate aging regime, the implications of our findings are important both for the understanding of aging in application and for the theoretical interpretation of the aging dynamics.By reference to the FD theorem, the consequence is that the properties governing the intermediate nonlinear physical aging of a system far from equilibrium are embedded in the thermal equilibrium fluctuations and can be predicted from these.This means that there is no fundamental difference between the intermediate nonlinear and the linear aging responses.Understanding physical aging is therefore intimately linked to characterizing and understanding the spectral shapes of linear responses and autocorrelation functions; a classical field where there has been important recent progress both experimentally [62,63] and theoretically [64].
The approach presented in this paper could prove useful for understanding the nonlinear response to electric fields.This is an active field [65][66][67] in which concepts from physical aging have been used successfully [66].
For future work it would also be interesting to see how far the description of physical aging in terms of linear response can be extended by including higher-order terms in the Taylor expansion of Eq. 8.This can hopefully lead to a complete picture of which samples and protocols exhibit aging governed by the same processes as those responsible of the linear alpha relaxation and which situations involve other processes and mechanisms [55][56][57].

MATERIALS AND METHODS
The study involves the glass-forming liquids VEC (99% purity) from Sigma-Aldrich for which data are shown in the main paper and NMEC (96% purity) from VWR for which data are shown in the Supplementary Materials.Both liquids were stored in a refrigerator at temperatures between 2 to 8 • C and used as received.
For each liquid a single sample was prepared for all the presented experiments.The sample cell was a plane-plate capacitor with a plate distance of 50 µm and a geometric capacitance of C geo = 16 pF.The cell was filled under ambient conditions and immediately mounted into a precooled cryostat.VEC was quenched to T cryo = 163 K, and NMEC was quenched to T cryo = 167 K, at which the samples were kept to equilibrate for a couple of days.The temperature of the main cryostat was constant at T cryo = 164 K for VEC and T cryo = 167 K for NMEC during the experiments that lasted almost 1 year for each sample.
The temperature control of the experiments was obtained by a microregulator integrated with the capacitor sample cell.The regulation was achieved by controlling a Peltier element in contact with a capacitor plate.Temperature was monitored with a negative temperature coefficient resistor placed inside one of the capacitor plates.A figure showing the sample cell with a microregulator can be found in the Supplementary Materials.Further details on the microregulator and the main cryostat are given in Ref. 37. The microregulator can change temperature by steps of a few millikelvin up to several Kelvin within seconds and keep the temperature constant with variations of less than 1 mK over weeks.All the temperature protocols shown, including the sinusoidal protocol, were achieved by making jumps in temperature with the microregulator.
The real and the imaginary part of the capacitance at 10 kHz was monitored during the entire experiment with a sampling rate of approximately one measurement per second.The measurements were performed using an AH2700A Andeen Hagerling ultra precision capacitance bridge.It is the combination of the fast and precise temperature control with the high precision of the bridge that makes it possible to measure aging in the linear limit.
The simulations used the Kob-Andersen 80/20 binary LJ mixture [60], which was simulated by means of standard N V T Nosé-Hoover dynamics [68] using the GPU-optimized software RUMD [69].A system of 8000 particles was simulated.In LJ units the time step was 0.0025.All pair potentials were cut and shifted at 2.5 times the length parameter σ ij of the relevant LJ pair potential (i, j = A, B).At the reference temperature T = 0.60 the potential-energy time-autocorrelation function was calculated as follows.First 10 7 time steps of simulations were carried out for equilibration.After that, the time-autocorrelation function was calculated using the Fast Fourier Transform.The temperature jump simulations were carried out by the following procedure applied for all starting temperatures.First, 5 × 10 5 time steps were spent on equilibration at the given starting temperature.After that a total of 5 × 10 8 time steps were spent on the production runs from which 1000 independent configurations were selected to serve as starting configurations for a temperature jump to T = 0.60.The Fig. 4  The temperature-control of the microregulator was calibrated with the temperature-specific microregulator voltage tracked during the measurements of spectra in the temperature range between 163 K and 180 K. Then the cryostat-controlled temperature was adjusted to T cryo =164 K and the micro-regulator was activated and set to T = T M R =164.5 K.When the sample environment (cryostat) is held at T cryo , the microregulator controls the exact sample temperature.Details on the microregulator setup are shown in Fig. S2.The microregulator can be tuned to ensure fast temperature equilibration, which is achieved after approximately 4 s for all temperature jumps performed on VEC.The microregulator ensures a long-time temperature stability better than 1 mK.
FIG. S2.Schematic drawing of the capacitor (dielectric cell) and microregulator mounted in the main cryostat [36].The liquid is inserted into the 50 µm gap between the disks of the dielectric cell (light blue).The Peltier element heats or cools the sample, depending on the direction of the electrical current powering the element.Details of the cryostat setup and the microregulator are given in Refs.37 and 70.

c. Overview of the loss-capacitance data
The dielectric response was measured with an Andeen-Hagerling ultra-precision capacitance bridge (model 2700A).The capacitance, C, and the dielectric loss, tan δ, were tracked at the frequency ν =10 Hz with a time resolution of about 1 s.
The data in the main paper was reduced in its resolution by averaging the data over linear and logarithmic bins, with exeption of the initial data points of individual data sets.This averaging does not affect to visual appearance of the data and does thus not affect their evolution or smoothness, but only serves for decreasing the size of resulting image files.The number of data points for individual-jump data is in the order of 100 after reduction of the resolution, decreasing the number of data points by up to three orders of magnitude.The most intense zooms in Figs. 2 and 3 in the main paper show data at full resolution.
Individual temperature jumps with amplitudes ranging from |∆T | =2 mK to 3.5 K were studied, as well as doublejump protocols with ∆T 1 =-100 mK, +100 mK, +1 K, +2 K or +7 K and ∆T 2 = −∆T 1 /2, and multi-jump temperature protocols resembling an overall sinusoidal temperature variation with amplitudes of 100 mK and 3 K.The full temperature protocol is shown in Fig. S3 together with the corresponding loss capacitance.The loss capacitance data are plotted with the original, as-measured resolution.
The main paper only considers the case in which the measured quantity X is the real part of the capacitance at 10kHz, so there we denoted the normalized relaxation function by simply R(t).The present document discusses also the case of X being the imaginary part of the capacitance for which reason we define the general normalized relaxation function as follows The single temperature jumps of the loss capacitance are plotted in Fig. S4 as functions of the logarithm of the time that has passed since the initiation of the jump, denoted by t − t i .To visualize the linearity of the jumps from Fig S4b, the normalized relaxation function is plotted in Fig. S5.The collapse of the data of all jumps with temperature amplitudes below 100 mK demonstrates a linear-response regime (at the experimental accuracy).As jumps of 2 mK show a considerably lower signal-to-noise ratio because the setup is here pushed closer to its limit in temperature control, only jumps with |∆T | ≥10 mK are depicted.Also the 50 mK down-jump is omitted in this plot.In Fig. S6, experimental data and predictions are plotted following double-jump temperature protocols in the linear regime.One of the realizations in Fig. S6a is included in the main manuscript; the accurate match with the second realization confirms the achieved results.The double-jumps in Fig. S6b follow an inverted Ritland-Kovacs crossover protocol, starting with a temperature up jump followed by a down jump.
Note that the predictions capture the details of the responses with a high accuracy: While the second jump is initiated at precisely the same point for the two realizations in Fig. S6a, a difference in the initiation of the second jump of the two realizations shown in Fig. S6b is observed.This leads to different levels in the short-time response at T = T 0 +50 mK, which are captured with a high accuracy by the predictions.
The predictions and experimental responses to single temperature jumps with amplitudes larger than 100 mK, i.e., nonlinear individual jumps, are plotted in Fig. S7.In addition to the data of the main manuscript, three additional jump amplitudes are presented here: |∆T | =0.5 K, |∆T | =1.5 K, and |∆T | =2.5 K.The details of the aging responses are captured accurately for all down jumps, and for up jumps with |∆T | ≤1 K, while the predictions show increasing deviations from the experimental data for larger up jumps.This signals a beginning deviation from the TN formalism that is known to work best for relatively small temperature variations.One possible scenario that qualitatively matches the observed behavior is the initiation of a heterogeneous growth process, as observed for very large temperature up-jumps for ultrastable vapor-deposited glasses.In Fig. S8, experimental data and predictions are plotted following double-jump temperature protocols in the nonlinear regime.In addition to the nonlinear double-jump that is included in the main manuscript, two more nonlinear double-jumps are presented here with temperature amplitudes ∆T 1 =+1 K in Fig. S8a and ∆T 1 =+7 K in Fig. S8c.Panels a to c show experimental data and predictions without any corrections.
From these three realizations of nonlinear double-jumps, only the 7 K double-jump exceeds the temperature regime for which equilibrium-data can be interpolated; its predictions are based on extrapolations far (3 K) below the known equilibrium values.Thus it is not surprising that the prediction for the 7 K double-jump deviates from the experimental data.To check under which parameter conditions predictions and experimental data coincide, we varied both the initial plateau of the second jump and the clock rate of the initial jump after making a minor adjustment of the overall data by shifting it by -6.5 fF along the y-axis to match to zero in the long-time limit of the second jump.The outcome is shown in Fig. S8d.The glassy contribution of the second jump was forced to match the data by multiplying the extrapolated parameter by 1.2 (visualized as grey diamond in Fig. S16).For the extrapolation of the clock rate the Avramov function was applied instead of the Vogel-Fulcher-Tammann (VFT) extrapolation, differing by a factor of 1.72 compared to the VFT-extrapolation (visualized as a black circle in Fig. S19).See Sec.A 2 for details for interand extrapolations of response-related parameters.

Loss contribution of dielectric capacitance
Figures S9 to S11 show comparisons of experimental data in terms of the loss contribution of the dielectric capacitance, C", and predictions derived from the loss contribution of the response to a linear 50 mK down-jump.C" (t,T) The predictions of the loss response to individual temperature jumps larger than 100 mK collapse convincingly with experimental data for down-jumps and up-jumps of amplitudes up to 1.5 K as shown in Fig. S9.Predictions of up-jumps with higher amplitudes appear retarded when compared to the experimental data.This behavior is qualitatively comparable to observations made on the storage capacitance, but is more pronounced for the loss data.The measured capacitance at constant temperature shows -even when the sample is thermally equilibrated -a small but steady decline (see Fig. S12).This "drift" is described and analyzed in detail in this section.
The data presented in the main manuscript and the previous section are corrected for drift by subtraction of a linear drift, which is either determined on the basis of the individual jumps ("individual drift correction") or sets of jumps ("loop-wise drift correction"), as described in the following.

Individual drift correction
For each jump that ends in equilibrium (individual jumps or second parts of double jumps), the region of constant slope in the long-time regime of the data is determined manually, typically involving of the order 10 4 data points.From a linear regression of these intervals according to the drift is identified as the slope r X (for the storage and loss response, separately).

Loop-wise drift correction
As for the individual drift correction, the region of constant slope in the long-time regime is determined manually for jumps that are assumed to have achieved equilibrium.For the loop-wise drift correction, subsequent jumps of same temperature amplitude |∆T | are regarded as a set of jumps, called a loop, typically consisting of four consecutive jumps (from or to T b by ±∆T ) or two consecutive jumps (from or to T b by |∆T |).The drift is determined from the slope s 1 of linear fits based on (for the response of the storage and loss component, separately).Thus we assume that all jumps of a loop may have different response "offsets" s 2 to s 5 , but share the same "drift" slope s 1 .
Based on the individual correction and the loop-wise correction, the raw time-dependent response data are corrected by subtracting the cummulated drift-contribution since the initialization of the experiment.The resulting corrected data sets are shown in Fig. S13a for the storage response and in Fig. S13b for the loss response.
In the below sections A 2 b to A 2 d it is shown that the difference in the data based on these two different driftcorrection methods is minor.This demonstrates the validity of both approaches.Since the equilibrium level at T ref is more stable for loop-wise corrected data, this correction has been applied to all linear data.Nonlinear data are corrected with the individual drift-correction because the high-temperature equilibrium values are more steady in this case.
The drift rate r X , determined from the individual drift correction, is plotted against the final temperature T b in Fig. S14.This figure shows temperature-dependent behavior for both the storage and the loss contribution; the evolution is more pronounced and systematic for the C -data.Here, the drift rate is fairly constant for T b ≤ T ref +1.5 K and increases with increasing temperature.For the loss data, r C" shows a weak increase towards positive rates where it becomes steady toward higher T b .The cause of this drift may be a temperature-dependent geometric change of the setup or some components thereof, an effect that is estimated to be weak compared to the thermally caused changes in the dielectric permittivity of the investigated materials.

b. Equilibrium values of the measured quantity Xeq
A sample that is held at constant temperature after a temperature jump is regarded as in equilibrium when the response levels off to a constant value at long times.In that stage, the physical aging process is concluded and the material reflects the properties of the supercooled liquid.In the raw experimental data, this long-time limit is characterized by a drift, as described in the previous section.The equilibrium values of the storage and loss response at a given temperature T b are determined from the drift-corrected data, defined as the mean value of the corrected data in the interval used to determine the drift.The equilibrium data were plotted against temperature and interpolated by a non-linear polynomial function as shown in Fig. S15.These data depend only weakly on the applied drift-correction method, as is evident from the parameters in table S1.The individually drift-corrected data and fits were chosen to represent the temperature dependence of X eq .On the basis of low-temperature data, a linear extrapolation was made that allows for estimating X eq for T ≤163.5 K.   S2.Parameters for the extrapolation along Ti and T b of individual and loop-wise corrected data for the glassy response contribution, R gl X , of the storage and the loss parts of the capacitance.

c. glassy contribution of the measured quantity
When plotting the normalized relaxation function, the level of the data just after the initiation of the jump, i.e. at small t − t i , is significantly below unity.This is due to the elastic contribution of the material and, possibly, also to secondary relaxation(s) occurring on time scales shorter than the time resolution of the experiment, i.e., at t − t i <4 s.This initial contribution is termed 'glassy contribution' in the following; it depends on both the initial temperature, T i , and the final temperature,T b , of a jump.The temperature dependence is fitted in the temperature ranges S2 for the parameters).Values connected to T i > T ref +1 K are excluded from the fit as the initial plateau is not sufficiently captured for jumps from elevated temperatures, resulting in an excess value of R gl X as it includes a relaxational contribution.The deviations between fits of individually corrected and loop-wise corrected data are small; thus the individually corrected data represents the temperature dependence of R gl X in Fig. S16.Avramov function: log 10 γ eq (T ) = w 1 − (w2/T ) w 3 ln (10) with w 1 =10.4,w 2 =275 K, and w 3 =6.89.
For the calculation of predictions, the value of γ eq was extrapolated by use of the VFT-function.For two individual jumps, A and B, the relation between their incremental time steps, dt * A and dt * B , is described by Eq. (s8).If one of these jumps ends at the reference temperature T ref , the normalized equilibrium clock rate can be written as γ eq,N (T b, B ) = γ eq (T b, B )/γ eq (T b, A = T ref ).Using the extrapolation of gamma values based on spectral data, we can check how well a fit of γ eq,N (T b, B ) relates to it.In this case the normalized equilibrium clock rate γ eq,N (T B ) is the only fitting parameter.The normalized equilibrium clock rate that is derived from the fit of storage data of individual jumps is plotted in Fig. S19 together with the spectra-based clock rates shown in Fig. S17, and which are normalized by the extrapolated value of γ eq at T ref .In order to calculate the prediction for a single temperature jump, the normalized response, R X (t), to a different single temperature jump is used as basis.This R X (t) is transformed back to its response amplitude in the duration of its time intervals to match the to-be-predicted data.This transformation utilizes the material-time aging formalism involving inter-and extrapolations of γ eq , Λ, and X eq , as well as the inter-or extrapolation of the glassy contribution R gl X .
The basis data set for all predictions of this work is a linear down-jump of 50 mK amplitude toward T ref .Note, however, that one is not restricted to linear jumps; the formalism allows for making predictions on the basis of any jump, linear or nonlinear.For VEC, only predictions for the storage data are presented.
In the following, details of the calculation of the various predictions are given.To clarify which parameters or data are connected to the linear data sets that serve as a basis for the predictions, the notation R lin is used for the normalized relaxation function for the linear jump of 50 mK amplitude that serves as a basis for the prediction, while R X is the normalized relaxation function for the to-be-predicted data set.As the basis data is a linear data set, its time can be referred to as the material time ξ.The individual data points of the data set are addressed by the index k.After the transformation of the material time to match the to-be-predicted data, the transformed time-data is referred to as the predicted time, t pred .The time data connected to the to-be-predicted data are referred to as laboratory time, t lab .

a. Predictions: Single temperature jumps
The following procedure describes the calculation of a prediction for an individual jump from T i to T b .The response amplitude of the prediction, ∆X(ξ k ), is calculated from the normalized relaxation function R lin (ξ k ) by with the overall response amplitude ∆X = X eq (T i ) − X eq (T b ).
The time interval of the prediction, dt k , is calculated according to the material-time aging formalism by with . The predicted time is given by b. Predictions: Double temperature jumps The prediction of the response to a temperature double jump can be understood as the prediction of an individual jump for the initial jump of the double-jump protocol, and a superposition of two individual jumps made in material time for the second jump of the protocol.For a visualization of this procedure, see Fig. S24.
The calculation of each of the double-jump predictions can be broken down into the two steps described above, i.e., the calculation of the response amplitude and the calculation of the time steps.Starting with the initial jump, the response amplitude of the prediction can be calculated as in Eq. (s11).Then, the time steps are transformed from material time to the time of the experiment as described in Eq. (s12).However, this is not done for the complete jump, but is interrupted at the point where the second jump is initiated.This is when the predicted time, t pred , is equal to the last point in laboratory time of the initial jump, t lab , of the to-be-predicted data set, just before the second jump is initialized.
Next the prediction for the second jump can be calculated.Since the point in time of the jump is known (t obs ), it is also known at which material time ξ the jump occurs, which is denoted by ξ 1 .As the prediction of the second jump is calculated from the superposition of the predictions for the individual jumps in material time, the response of the total prediction is given by The transformation from material time to t pred is analogous to the approach described for the individual jumps in Eq. (s12), where ∆X is replaced by ∆X tot from Eq. (s14).

c. Predictions: Multiple temperature jumps
The prediction of a response that involves multiple jumps is very similar to the prediction of a double jump: the approach iteratively determines the response amplitude of a prediction as a superposition of the individual responses in material time.This is again followed by the transformation of the time step from material time to t pred in order to match the laboratory time of the to-be-predicted experiment.The only difference is that the superposition in material time involves several individual predictions: b. Storage and loss capacitance data for NMEC After the measurements of the above spectra, the microregulator was activated and set to T = T ref =166.1 K.By means of the Andeen-Hagerling capacitance bridge the capacitance, C, and the dielectric loss, tan δ, were tracked at the frequency ν =10 Hz.Individual temperature jumps with amplitudes from |∆T | =5 mK up to 3 K were initiated, as well as double jumps with ∆T 1 =+50 mK and ∆T 2 = −∆T 1 /2.An overview of the temperature protocol is shown in Fig. S21 together with the corresponding storage and loss capacitance.The first four jumps were measured with the dielectric setup used for the measurements of the spectra in Fig. S20.These data emphasize the significantly higher accuracy of measurements with the Andeen-Hagerling AH2700 used for obtaining the aging data.After approximately 35 weeks of consecutive measurements, a failure in the cooling unit of the the cryostat unfortunately terminated the measurement sequence.
The analysis of the capacitance data for NMEC mainly follows the procedure described for VEC in the previous section.However, for NMEC only a loop-wise drift correction was performed, individually corrected data was not generated.
Note the difference in signal-to-noise ratio for the conventional dielectric setup and the Andeen-Hagerling instrument that becomes obvious in the fourth jump of the data set (T ref to T ref +3 K, red color), where the measurement setups were switched.The response to individual temperature jumps is plotted against logarithmic time in Fig. S22.The initial 10 s of each jump are clearly still influenced by the temperature regulation and were excluded from the analysis.After this thermalization, the initial plateau is only captured for jumps from T < T ref to higher temperatures (purple colors in Fig. S22a); for all other jumps a significant amount of relaxation has taken place, implying that the data set on NMEC does not allow for an analysis of the glassy contribution.Thus, this analysis was omitted in case of NMEC.
Fig. S22b depicts the individual temperature jumps with temperature amplitudes between 5 mK and 100 mK that correspond to the data magnified in the insets of Fig. S21.In Fig. S23 the linear aging data are plotted in the normalized representation.Data with amplitudes down to 10 mK collapse as expected for linear relaxation.The responses to 5 mK-jumps show deviations from the master curve, even though the data still follow the evolution of the collapsed data.This deviation is clearly influenced by thermal fluctuation.In order to predict the response to temperature jumps larger than 100 mK not only the equilibrium reponse value X eq and the response to a linear jump (50 mK up-jump to T ref ) is used, but also the single parameter Λ and the equilibrium clock rate γ eq .As for VEC and as used for the prediciton of linear double jumps, X eq was determined for the storage and loss capacitance analogous to the procedure described for VEC in section A 2 b.
The parameter Λ for NMEC was determined by the integral criterion as described in section A 3 b for VEC.Fig. S25 visualizes the minimization of the sum of integrals I 1 + I 2 for the case of storage data, resulting in a value of Λ C = 6.80 pF −1 .The equilibrium clock rate is based on an extrapolation from the fit of relaxation times derived from dielectric spectra to a Vogel-Fulcher-Tammann function, analogous to the description given in Sec.A 2 d for VEC.However, clock rates determined from fits to the storage response of individual jumps at the highest investigated temperatures do not collapse with the VFT-extrapolation based on clock rates derived from spectra, as visualized in Fig. S26.
The storage contribution of the capacitance for individual jumps with amplitudes larger than 100 mK, i.e., nonlinear individual jumps, are plotted in Fig. S27 together with predictions.As described before, predictions are based on a linear response data set R lin , equilibrium response values, C eq , the equilibrium clock rate, and the single parameter, Λ C .Again, predictions are based on a 50 mK up-jump as used for the linear case as in Fig. S24.Equilibrium response values are interpolated analogous to the procedure described for VEC, and also the single parameter Λ is derived in analogy as shown in Fig. S25.In order to take varying levels of the glassy response contribution into account, the as-measured levels were applied instead of making an interpolation, as the initial plateau value at short times is not well resolved for most individual jumps for NMEC.As shown in Fig. S26, the values of the clock rate that are derived from fits of the individual jumps do not follow the tested extrapolations as nicely as the for VEC.In Fig. S27 the outcome for predictions based on fitted (a) and extrapolated clock rates (b) are visualized.While the fitted clock rates yield a better collapse between predictions and experimental data, the extrapolation of clock rates by the VFT-function yields qualitatively good results but more deviations than VEC-predictions, with strongest differences for jumps from the reference temperature upwards, similar to the trend observed for VEC.

FIG. 1 .
FIG.1.Overview of the temperature protocol and the raw data of the full experiment on VEC.A The experimental protocol realized by temperature modulations around the reference temperature 164.6 K (upper panel) and the real part of the measured capacitance C (10 kHz) (lower panel), both plotted as functions of time on a linear scale.Jumps larger than 100 mK are colored while jumps of 100 mK or less are depicted on a grey scale; a selection of the latter is shown in the inset.The sinusoidal temperature modulations that are also studied (see below) are not resolved in this figure where they appear as turquoise thick vertical lines.B The capacitance C (10 kHz) data plotted as functions of the logarithm of the time t that has passed after each jump.C Magnification on the jumps of magnitude 100 mK or less.

FIG. 2 .
FIG.2.Data for the real part of the capacitance, C (10 kHz), from small-amplitude temperature modulation experiments on VEC along with predictions demonstrating that the response is linear.A Normalized relaxation function (Eq. 3) of single temperature jumps of amplitude 10 mK-100 mK around the reference temperature T0 = 164.6K.All data collapse as predicted for linear relaxation.B Data from a small-amplitude temperature double jump starting at T0 = 164.6K and jumping first by -100 mK and then by +50 mK (colored curves) along with the prediction according to Eq. 4 (black dashed lines).The upper panel shows the full experiment on a linear time scale, the lower panel shows the data on a logarithmic time scale that sets the time of the beginning of each temperature jump to zero.The inset illustrates the temperature protocol.C Sinusoidal small-amplitude temperature protocol (red points) and data (turquoise points behind the black points).The prediction based on the linear-response formalism (Eq.5) is shown as black points.The deviations between prediction and data that can barely be discerned in the main figure can be seen in the magnification.
FIG. S1.Dielectric spectra for VEC.Storage (a) and loss (b) dielectric permittivity, ε and ε", as functions of frequency measured by a custom-built frequency generator and a commercial LCR meter in the temperature range 163 K-180 K.The inset in panel a depicts the skeletal formula of VEC.
FIG. S3.Overview of the temperature protocol of consecutive physical-aging experiments realized by temperature jumps facilitated by an NTC-thermistor-regulated Peltier element around T ref =164.6 K, and of the loss capacitance, C (ν = 10 kHz), as functions of time (in seconds on the lower and weeks on the upper x-axis).Nonlinear aging data are colored while linear data are depicted on grey scales.The insets show single jumps with temperature amplitudes of 100 mK or less with a zoom on the y-axes.
FIG.S8.Experimental data and predictions of nonlinear double jumps on VEC involving (a) temperature jumps of ∆T1 =+1 K and ∆T2 =-0.5 K, (b) ∆T1 =+2 K and ∆T2 =-1 K, (c) ∆T1 =+7 K and ∆T2 =-3.5 K for two realizations.The predictions are based on inter-and extrapolation of the glassy contribution X gl , the equilibrium response Xeq, Λ, and the equilibrium clock rate γeq.Note that for the predictions of the double jump with ∆T1 =+7 K shown in panel (c) the aforementioned quantities had to be extrapolated, while the predictions of double jumps in panels (a) and (b) are based on interpolations.For more details on how the predictions are derived, see section A 4 b.No corrections were applied to yield the predictions shown in panels (a) to (c).In panel (d), the data were shifted by -6.5 fF along the y-axis to match to zero in the long-time limit of the second jump, and corrections for the clock rate and the glassy contribution of the second jump were applied.

2 .
FIG. S10.Experimental data and predictions of linear double jumps on VEC involving (a) temperature jumps of ∆T1 =+100 mK and ∆T2 =-50 mK and (b) ∆T1 =-100 mK and ∆T2 =+50 mK, each showing data for two realizations.Upper panels show data plotted against linear time, middle and lower panels depict data as functions of the logarithm of the time that has elapse since the initiation of a jump at t = ti.No corrections were applied to yield the match between experimental data and prediction.

9 FIG
FIG. S12.Overview of the consecutive physical-aging experiments based on the temperature protocol presented in the main manuscript without drift correction for (a) the storage and (b) the loss contribution to the measured capacitance.
FIG. S15.Equilibrium responsesXeq in form of equilibrium values of (a) the storage and (b) the loss capacitance as functions of temperature for individually-corrected data.Full lines correspond to nonlinear interpolations of the data; dashed lines represent linear extrapolations based on the data between 163.5 K and 164.5 K.
FIG. S16.The glassy response contribution, R gl X , for jumps starting at T ref plotted as a function of final temperature, T b (panel ((a)), and for jumps ending at T ref plotted as a function of initial temperature, Ti (panel (b)).The full lines represent the interval used for the fit, while the colored dashes lines represent the extrapolation over the full temperature range.The gray dashed line in the lower panel of (b) after a temperature jump represents the extrapolation of the glassy response contribution for jumps from temperatures Ti to Tb =164.1 K, i.e., R gl X (Ti, T b = T ref −0.5 K).The gray diamond reflects the value of R gl X that is applied for the second jump of the 7 K double-jump protocol in Fig.S8d.All values shown in the lower part of panel (b) can be fitted by a second-order polynomial (not shown).
FIG. S17.Temperature dependence of the inverse time scale derived from spectral loss-peak positions, γeq, and extrapolations based on the Vogel-Fulcher-Tammann (VFT), parabolic, and Avramov functions.

1 ) 4 .
FIG. S19.Logarithmic normalized equilibrium clock rate, γeq,N as a function of the temperature.Dashed lines correspond to extrapolations based on the clock rates determined from spectral data (grey cirles), crosses correspond to clock rates from fitting the storage contribution of individual jumps.The inset shows a zoom of γeq,N derived from the individual jumps and extrapolations from the spectral data.

FIG. S21 .
FIG. S21.Overview of the temperature protocol and raw data of the full experiment on NMEC after a loop-wise drift correction.The temperature protocol of experiments realized by modulations around T ref = 166.1 K is plotted in the upper panel together with the real and loss part of the measured capacitance C(ν =10 Hz) in the lower panels, each as functions of time on a linear scale.Jumps larger than 100 mK are colored while jumps of magnitude 100 mK or less are depicted on a gray scale.The insets show details of jumps with |∆T | ≤100 mK and share the same time-axis as the main panels.
FIG. S22.The storage and loss capacitance data, C (ν =10 kHz) and C"(ν =10 kHz), plotted as functions of the logarithm of the time that has passed after the initiation of each jump at ti. Panel (a) depicts jumps with |∆T | >100 mK and panel (b) is a zoom on jumps of magnitude 100 mK or less.
FIG.S23.Normalized relaxation function for the real and loss parts of the capacitance data, R C and R C" , of single temperature jumps with amplitudes between 5 mK and 100 mK.
Fig. S24 visualizes the superposition of individually predicted jumps to yield the overall prediction for the doublejump temperature protocol.This example covers a double jump in the linear regime with ∆T 1 =+50 mK and ∆T 2 =-25 mK.In (a) the storage and loss capacitance are plotted against linear time.(b) reflects distinct points in the thermal protocol, such as the initial equilibrium, the glassy response contribution, the initiation of the second jump, and the final equilibrium, in a plot of the dielectric response as a function of temperature.Storage data from (a) are presented in the lower panel of (c) in comparison to the calculated prediction.(d) shows a comparison of storage and loss capacitance from (a) to predictions on a logarithmic time scale.More details are given in the figure caption.a) b) FIG.S25.Determination of Λ.The sum of integrals, I1 + I2, plotted as a function of Λ derived from loop-wise corrected data of a pair of temperature jumps of 100 mK amplitude towards the reference temperature, T ref .

FIG. S27 .
FIG. S27.Experimental data of nonlinear individual jumps on NMEC involving temperature jumps between 1 K and 3 K.The predictions (black lines) are based on different parameter sets: a) predictions are based on the glassy contribution as observed for the experimental data for individual jumps and with clock rates from fits to the experimental data for individual jumps.b) predictions are based onthe glassy contribution as observed for the experimental data for individual jumps as in a), but clock rates are extrapolation by the Vogel-Fulcher-Tammann function from loss-peak positions determined from spectral data.

TABLE S1 .
corr.loop-wise corr.indvl.corr.loop-wise corr.Parameters for nonlinear interpolation and linear extrapolation of individual and loop-wise corrected equilibrium data, Xeq, of the storage and the loss parts of the capacitance.
FIG. S18.Determination of Λ for storage and loss contribution of the data.The sum of integrals, I1 +I2, plotted as a function of Λ derived from individually corrected data of a pair of temperature jumps of 1 K amplitude towards the reference temperature, T ref .
c. Check on the equilibrium clock rate