nature physics Article https://doi.org/10.1038/s41567-023-02125-0 Real-time microscopy of the relaxation of a glass Received: 13 January 2023 Marta Ruiz-Ruiz1,2, Ana Vila-Costa 1,2, Tapas Bar 2, Cristian Rodríguez-Tinoco 1,2 , Marta Gonzalez-Silveira , 1,2 Accepted: 4 June 2023 Jose Antonio Plaza 3, Jorge Alcalá 4, Jordi Fraxedas 2 & Published online: 13 July 2023 Javier Rodriguez-Viejo 1,2 Check for updates The understanding of the dynamics of a glass above its devitrification temperature remains incomplete. Here, we build a spatio-temporal map of the relaxation dynamics of a highly stable glass into its supercooled liquid using real-time atomic force microscopy imaging. This methodology enables direct visualization of the progression of the liquid phase and clarifies and quantifies the presence of localized fast mobility regions separated by giant length scales. Our data establish a clear correlation between dynamic length and time scales in glasses. This approach may also be applicable to unveil the microscopic structure and dynamics of other glass-forming systems with much shorter length and time scales, including liquid-cooled glasses. Glasses are non-equilibrium materials arrested during cooling in met- though the associated length scale has not been directly measured by astable configurations. The temperature that marks the transition experiments yet15–17. Vapour-deposited stable glasses appear as model from the ergodic state, the liquid, to the non-ergodic one, the glass, is systems to explore this scenario18–21 and recent simulations point in known as the glass transition temperature, Tg. A hallmark in liquid and this direction22. glassy dynamics is the recognition that dynamic heterogeneities are Since direct spatial visualization of equilibrated regions above at the core of the slowing down of the dynamics and are responsible Ton is extremely challenging due to the subtle structural changes over for the glass transition and its temperature dependence1,2. Dynamic very small distances, we follow a different strategy based on the local heterogeneities are characterized by clusters of regions of atoms or mechanical instabilities that the liquid regions generate on a rigid molecules with correlated mobility that grow as the temperature is ultrathin layer grown on top of the glass. Well-defined surface undu- decreased. The existence of these clusters of mobility has been inferred lations can be produced in several ways23,24, leading to self-organized from previous experiments1,3–11, but it is experimentally challenging to wrinkling patterns with a rich variety of morphologies depending on identify them or their length scales, although it is recognized that the the particularities of the system25,26. One approach to induce surface broad distribution of time scales in a glass can be associated with the wrinkling is by capping a soft thin-film material (typically a polymer, presence of mobile and immobile regions12,13. It is generally accepted but also a small-molecule organic glass) spin-coated on a rigid sub- that the relaxation of the glass into the supercooled liquid (SCL) at the strate by a thinner film of a metal layer24 or by another organic layer devitrification temperature on heating, Ton, happens through a gradual with a higher Tg (ref. 27). The instability can be initiated by annealing softening of the glass across its entire volume with a correlation length the organic material into the rubbery state where wrinkling across of a few nanometres. However, we have recently shown this view is the whole surface appears due to the development of compressive not unique and can be changed by accessing a temperature and/or stresses induced by differences in thermal expansion coefficients28,29. time regime above Ton where glasses transform into the liquid by the Wrinkling can be locally induced at specific sites of the film by localized formation of localized regions of liquid within a glassy matrix14, even surface modifications using focused ion beams30, swelling-induced 1 Departament de Fisica, Facultat de Ciencies, Universitat Autonoma de Barcelona, Bellaterra, Barcelona, Spain. 2Catalan Institute of Nanoscience and Nanotechnology (ICN2), CSIC and BIST, Campus UAB, Bellaterra, Barcelona, Spain. 3Instituto de Microelectrónica de Barcelona, IMB-CNM (CSIC), Campus UAB, Bellaterra, Barcelona, Spain. 4InSup, ETSEIB, Polytechnic University of Catalonia, Barcelona, Spain. e-mail: cristian.rodriguez@uab.cat; javier.rodriguez@uab.cat Nature Physics | Volume 19 | October 2023 | 1509–1515 1509 Article https://doi.org/10.1038/s41567-023-02125-0 a b Temperature Isotherm Tg + 16 K TCTA–13 nm Qh = 5 K min–1 Tg TPD–63 nm Tg – 10 K TCTA–13 nm Si substrate Troom Q = 5 K min–1 Time c 125 min 170 min 200 min End iii 2.0 20 d 15 Intensity (a.u.) 1.5 Height (nm) 10 10 1.0 5 0.5 0 0 d = 912 ± 70 nm 0 –5 –10 0.2 0.3 0.4 0.5 0.6 0.5 1.0 1.5 0.5 1.0 1.5 2.0 2.5 3.0 3.5 –4.0 –3.5 –3.0 –2.5 –2.0 Length (µm) Length (µm) Length (µm) k (nm–1) Fig. 1 | Representative images of the highly heterogeneous transformation of Scale bars, 4.0 μm. d, Power spectral density (PSD) distribution of the fully a stable TPD glass. a,b, Sketch of the trilayer structure (a) and of the annealing transformed sample indicating the prominent wavelength of the wrinkled protocol (b) in the AFM setup. c, AFM snapshots obtained at 349 K after time distribution. PSD has been calculated with WxSM software with no associated t showing the formation and progression of local wrinkles induced by the uncertainty. The value of d is determined from the maximum of PSD with an error transformation of the middle TPD layer into a SCL. The bottom graphs show the corresponding to the separation between two adjacent points. line profiles of the selected regions. A white arrow marks the pattern considered. stress through toluene absorption31 or local heating with an external is shown in Fig. 1a,b. The TPD layer is vapour-deposited at a sub- source, such as a laser32. strate temperature of 0.85Tg1 (285 K) and has a fictive temperature37, Here, we take advantage of a previously established proto- Tf = 292 ± 2 K, 40 K below its Tg1. Details of the growth of the organic col to induce bulk ‘melting’ in thin-film glasses33–35 by capping the layers and the experiments are given in the Methods. After growth, the organic film (N,N′-bis(3-methylphenyl)-N,N′-diphenylbenzidinem, thin-film stack (13 nm TCTA/63 nm TPD/13 nm TCTA/silicon, Fig. 1a) TPD, Tg1 = 333 K measured at 10 K min−1) with two ultrathin layers of an is taken to an AFM or optical microscope equipped with a heating organic glass with a higher Tg (tris(4-carbazoyl-9-ylphenyl)amine or stage. Most results involve AFM due to its superior spatial resolution. TCTA, Tg2 = 428 K, 10 K min−1). The mechanical instabilities produced In AFM, the temperature is raised to a working value of 349 K, Tg1 + 16 K, during isothermal treatments above Tg1 are measured in situ using in sequential steps as shown in Fig. 1b to allow for equilibration of the atomic force microscopy (AFM) and optical microscopy. We address the AFM setup (Methods). At 349 ± 1 K, the transformation time of TPD into continuum mechanics finite element modelling (FEM) of the wrinkling the SCL is around 270 min (ttrans ≅ 2 × 106τα, where τα is the relaxation phenomenon to demonstrate that the local surface wrinkles appearing time of the equilibrated liquid at that temperature)38, taking t = 0 s at Tann > Tg1 (TPD) are induced by the equilibration of the intermediate as the time at which this temperature is reached. This time interval TPD layer at localized spots. The observed patterns are compatible with offers a convenient scale to follow the kinetics of the transformation a time-dependent number of initiation sites that appear to propagate in real time, given that each AFM scan takes around 4 min to complete radially. With these data, we construct a spatio-temporal, microscopic (Methods). Figure 1c shows some representative in situ AFM images map of the relaxation of a molecular glass, delivering a clear picture of obtained at T = Tg1 + 16 K in the same spatial location for different times. the heterogeneous devitrification dynamics of glasses with enhanced The surface of the trilayer was initially very smooth (RMS roughness stability. Imaging by AFM or optical microscopy provides a direct below 0.3 nm in 1 × 1 μm2) and did not show any surface undulation visualization of giant lengths scales between fast mobility regions and before reaching the annealing temperature. Figure 1c bottom repre- the growth of the liquid by dynamic facilitation. This methodology sents the apparent height of the surface topography across the line also offers an original approach to study the glass transition from a scans marked in the AFM images (Methods). The images reveal the localized and new perspective and bridge the gap with the microscopic emergence of surface corrugation that initially appears in the form of information accessible to computer simulations13,22,36. Gaussian-type protuberances with aspect ratios around 100:1, that is, A sketch of the trilayer organic structure grown on a 500-μm-thick a few nanometres in height and a few hundred nanometres in width. Si substrate with its native oxide (surface root mean square (RMS) The corrugation develops over time so as to form localized wrinkle pat- roughness, 0.13 nm in 1 × 1 μm2) together with the thermal protocol terns that grow radially from the initial seed. A set of AFM and optical Nature Physics | Volume 19 | October 2023 | 1509–1515 1510 Article https://doi.org/10.1038/s41567-023-02125-0 a b c 1.0 FEM AFM 0.8 Amplitude (a.u.) 0.6 0.4 0.2 0 0 200 400 600 800 1,000 x (nm) 0 1 d TCTA TPD TCTA Si –0.5 1.0 Fig. 2 | Surface undulation patterns and their modelling via finite element sectional view of the surface undulation region. The neo-Hookean model with simulations. a, Schematic showing the growth of an SCL region (in green) material parameters C = 3.71 × 106 Pa and d = 5.58 × 10−8 Pa−1 is assumed in the from a nucleation site within the TPD interlayer towards the top and bottom simulations. c, Comparison between the simulation results from finite element TCTA layers (growth direction marked with green arrows). This results in the modelling (FEM) and the experimental AFM measurements concerning the development of the first surface undulation triggered under the applied thermal shape of the first undulation during the early propagation stages of the liquified stresses (marked with black arrows). Notice the periodic onset of further front. AFM data is for an emerging liquid nucleus, with no associated uncertainty secondary undulations as the cylindrical-shaped SCL region extends radially values. d, Comparison between a simulated wrinkled pattern of θ ≅ 1,000 nm and mechanical instabilities develop. b, Simulated normalized out-of-plane under the assumption of the neo-Hookean model with material parameters displacement (top) isotropic view of the whole simulated structure, (bottom, C = 3.71 × 106 Pa and d = 5.58 × 10−8 Pa−1, and the AFM image of a typical pattern in left) isotropic view of a zoomed region containing the first surface undulation 13/63/13 nm trilayers. Scale bars, 1.2 μm. of diameter θ ≅ 250 nm within the TPD interlayer and (bottom, right) a cross- microscopy images of the transformation of a 13/63/13 nm trilayer localized patterns evolving with time. TPD samples grown directly on are shown as videos in Supplementary Videos 1 and 2. We note that weakly interacting native oxide–Si substrates and capped only on top the onset of the local perturbations could not be identified by AFM, by TCTA do show wrinkling at the whole surface from the beginning so these measurements are not sensitive to the initial formation of the of the experiment due to the transformation of TPD into the SCL by a liquid droplet within the film and only monitor its appearance when propagation front that starts at the TPD–SiO2 interface (Supplemen- the liquid region is larger than the thickness of the middle layer (sche- tary Fig. 2). This corrugation behaviour is consistent with the thermal matics in Fig. 2a). This limits the spatial resolution to the thickness of wrinkling commonly observed in metal and polymer bilayers39. this layer. The peak-to-peak amplitude of the central part evolves from 2 to 3 nm of the initial seeds to around roughly 30 nm (peak-to-valley) Mechanical analysis of wrinkle formation at the end (Fig. 1c), although the precise dimensions depend slightly Continuum mechanics finite element simulations are used to correlate on the selected wrinkled pattern. The dominant wavelength of the the formation of surface undulations with the existence of SCL regions fully transformed sample, obtained from the power spectrum analy- within the intermediate TPD layer. A distinctive feature of the investi- sis of the fast Fourier transform image, Fig. 1d, amounts to 912 nm gated mechanical instabilities compared with most previous thermal (Methods). The AFM images in Fig. 1 and Supplementary Videos 1 and 2 wrinkling studies is their local character, where the softened (liquid) show two main features: the appearance of new local corrugations sites slowly populate the interlayer. While the three-layer thin-film of the trilayer as time evolves and their radial growth. Corrugation is system thus behaves in accordance with linear thermoelasticity, cylin- irreversible and maintains its structure over months at room tempera- drical SCL regions arise in the middle TPD layer after some time above ture. The observed patterns are not associated with crystallization, Tg1. Surface undulation patterns are then induced due to the transfer grazing-incidence X-ray diffraction as well as transmission electron of the compressive thermal loads from the TPD interlayer to the stiffer microscopy imaging of the fully corrugated sample after the thermal and much thinner TCTA top layer, which remains in its glassy state. This treatment show that the trilayer maintains its amorphous structure is a mechanically unstable phenomenon triggered under the roughly (Supplementary Fig. 1). Capping the TPD middle layer at both sides 33-fold higher elastic modulus and roughly 88-fold smaller thermal is required to avoid surface induced transformation and to observe expansion coefficient of the Si substrate compared to the organic Nature Physics | Volume 19 | October 2023 | 1509–1515 1511 Article https://doi.org/10.1038/s41567-023-02125-0 layers. These surface undulation processes are schematically illustrated or (2) it is directly related to the intrinsic shape of the distribution in Fig. 2a, resulting in the experimentally observed wrinkled patterns. of relaxation times of the UG at this temperature, or to some mixed Concentrically ridged undulations are also produced as a function of effect. In Fig. 3b we plot a two-dimensional spatial map of the sample the relative layer thicknesses (Supplementary Fig. 3). with the position of the equilibration spots as they appear. The colour The onset of the first surface undulation (single ridge configu- scale marks the time from the appearance of the first liquid regions. ration) is concomitant to the vertical expansion of the cylindrical An initial period of about 2 h without any noticeable surface fea- SCL regions due to density (ρ) differences between the liquid and ture can be clearly discerned in Fig. 3a. This long transient is related to the ultrastable glass (UG) (ρ(UG) 2.8% higher than ρ(SCL) at Tg1 + 16 K the high stability of the glass, that is, some regions attain liquid mobility (refs. 40,41) and the compression induced by the surrounding TPD within this long time. To some extent, the large induction period could glassy interlayer. The capacity of the cylindrical SCL regions to with- also be associated to a slower growth velocity of the minuscule initial stand shear stresses at the outermost regions close to the top TCTA seeds localized inside the TPD matrix before they reach the TCTA cap- layer is then regarded as an essential element of the attainment of the ping layers due to pressure effects originating from density differences multiple periodic undulations associated with the mechanical instabil- between the UG and the SCL above Tg1, as reported recently in a numeri- ity of the thin-film system. The sudden inception of fully wrinkled sur- cal simulation study22. This effect may have a minor influence here faces with the partly capped, two-layer system (Supplementary Fig. 2) because of the small layer thickness but could become important for is indicative in that this capacity is maintained across larger surface thicker films. In the time interval between 0.9 and 1.25 × 104 s, the evolu- areas. Some viscous relaxation of the shear stresses is finally seen to tion of the number of liquid seeds (Fig. 3a) can be roughly approxi- occur at Tg1 + 16 K, tentatively reducing the heights of the undulations. mated by a linear function (dashed blue line in Fig. 3a) giving an average In our finite element simulations, we consider that the softened ‘formation frequency’ ν ≅ 6 ± 2 × 107 nuclei m−2 s−1. This magnitude TPD interlayer is a dense liquid that exhibits solid-like properties below translates into a mean separation distance of liquid droplets appearing a certain elastic limit. As we are essentially interested in assessing the in 1 s per unit area of ≅130 μm, using ξ(T) = (ν (1s))−1/d , where d is the onset of the surface undulations under sustained shear, it is assumed dimensionality of the system. Experimentally, the giant length scale that the cylindrical SCL regions deform though linear elasticity with of the bulk transformation was previously indirectly inferred from the discretionary levels of compressibility and small elastic modulus. lower bound of the cross-over length between front propagation and The softened interlayer is alternatively assumed to deform under bulk melting, as well as from calorimetric measurements that identified the non-linear neo-Hookean hyperelastic model, so as to capture a the appearance of isolated liquid regions15,16. Length scales of the order possible decrease in the initial elastic modulus at small strain levels of 1–4 μm have been suggested. To compare with previous estimates (Methods). Comparisons between the experimental and computational of the cross-over length in thin films, it is necessary to account for the results are shown in Fig. 2b–d. A cross-section of a scanning electron time dependence. That is, in 3,000 s the average distance between micrograph showing the surface undulation of a trilayer is shown in localized liquid regions roughly equals 2 μm (under the assumed Supplementary Fig. 8b. ν = 6 × 107 nuclei m−2 s−1 (Fig. 3a)), which falls within the same order of magnitude as the roughly 1.2 μm distance travelled by the front at Equilibration of liquid regions: giant length Tg1 + 16 K in the same time scale (considering vfront ≅ 0.4 nm s−1, as dis- scales cussed below). This confirms the length scale from the previous The images of Fig. 1c and Supplementary Videos 1 and 2 offer a detailed cross-over experiments15,16. view of the evolution of the glass following a temperature jump above Tg1, in which localized regions within the TPD glass having higher mobili- Relaxation map of the glass ties evolve towards the liquid state: a process referred to as softening We also follow the time evolution of the size of the liquid regions both and equilibration, subsequently consuming the glassy environment by from AFM and optical microscopy images (Fig. 4 and Supplementary dynamic facilitation. These images clarify previous experiments and Fig. 5). The methodology used to extract the radius of the liquid theories about the bulk transformation of stable glasses18–20,42. These domains and the transformed fraction of liquid from the AFM and studies, through indirect measurements, suggested that bulk melting on optical microscopy images is described in detail in the text accompany- these glasses proceeds through a ‘nucleation and growth’-type process ing Supplementary Figs. 12 and 13. A strong indicator that the wrinkles with giant length scales between the initial patches of liquid. Here, we originate from the presence of liquid in the TPD layer is found in the assume the first regions equilibrating into the liquid correspond to those radial propagation velocity of the instability region. We have analysed with the highest mobility and therefore we interpret the emergence many individual spots, such as those shown in Fig. 4a, rendering an of newly equilibrated regions at different times as clear evidence of average radial propagation speed at Tg1 + (16 ± 1) K of 0.3 ± 0.1 nm s−1. dynamic heterogeneity. Our results confirm that dynamic heterogenei- This is consistent with the growth front velocity of the SCL TPD at this ties in UGs are characterized by giant length scales of many micrometres temperature, v = 0.4 ± 0.1 nm s−1, evaluated independently through between fast mobility regions. Besides these localized spots, we also extrapolation of high and low-temperature calorimetric data (Supple- have indirect evidence that the entire glass matrix slowly softens during mentary Fig. 6). Thus, the wrinkle patterns propagate laterally at the the annealing above Tg1, but to the most part without reaching equili- speed of the liquid front, confirming a previous hypothesis17. Although bration. Calorimetric measurements show the change in the position the liquid signature on the mechanical instability looks radial on aver- of the endothermic overshoot of the highly stable glass matrix during age, we cannot infer whether mobility locally spreads isotropically or the annealing treatment (ref. 17 and Supplementary Fig. 4). through string-like, fractal, dynamics. From the wrinkled area, we can In Fig. 3, we represent the time evolution of the heterogeneous, estimate the liquid transformed fraction (Supplementary Fig. 13) that localized softening obtained from AFM. Figure 3a plots the cumulative follows a sharp sigmoidal shape, as shown in Fig. 4b, most probably due number of regions that enter the regime of equilibration and become to the continuous formation of new liquid seeds and its acceleration as liquid-like per unit area (black symbols). The inset is a plot of the newly time evolves. Figure 4c,d represents spatio-temporal maps of the evolu- developed liquid regions per unit area (black circles), and normalized tion of the liquid including its progression by dynamic facilitation taken by the area of the untransformed fraction (blue symbols) as a function at transformed liquid fractions X = 34 and 90%, respectively, and can be of time (Supplementary Information). This representation highlights interpreted as relaxation maps of the glass. As time progresses, regions an acceleration of the number of local regions entering the liquid of fast mobility spatially separated by giant length scales equilibrate regime over time. We foresee two possible scenarios to explain this towards the liquid and they ignite adjacent regions that also become behaviour: (1) acceleration originates from the softening of the matrix, more mobile by a facilitation mechanism. A direct comparison between Nature Physics | Volume 19 | October 2023 | 1509–1515 1512 Article https://doi.org/10.1038/s41567-023-02125-0 a b ×104 8 Number of liquid seeds (×1011 m–2) 1.4 3 ∆N (×1010 m–2) 6 1.3 4 1.2 Time (s) 2 2 1.1 0 5 10 15 1.0 1 Time (×103 s) 0.9 0.8 0 0 2.5 5.0 7.5 10.0 12.5 15.0 Time (×103 s) Fig. 3 | Formation of liquid seeds in the TPD interlayer during annealing in the inset is propagated from the error in the calculation of the transformed above Tg. a, Evolution of the cumulative number of nuclei per unit area. No fraction (Fig. 4b). b, A spatial map at 90% of transformation obtained from AFM uncertainty is considered for the absolute discrete number of nuclei. The inset images showing mobility regions where glass relaxes directly into the liquid. shows a number of new nuclei that appear in each AFM image per unit area (black The colour code reflects the time of formation; t = 0 marks the time at which the symbols) and per unit area of untransformed glass (blue symbols). The error bars temperature reaches 349 K. The first soft spots appear after 2 h. a 1.5 b Nuclei-1 100 Nuclei-2 Nuclei-3 Nuclei-4 80 Nuclei-5 1.0 Nuclei-6 Nuclei-7 Radius (µm) 60 Linear fit x (%) 40 0.5 20 0 0 3 3 3 3 3 4 4 0 1 × 10 2 × 10 3 × 10 4 × 10 5.0 × 10 1.0 × 10 1.5 × 10 Growing time (s) Time (s) c ×104 d ×104 6 1.2 6 1.5 1.4 1.1 7 7 1.3 1.2 1.0 Time (s) Time (s) 4 4 1.1 2 5 0.9 2 5 1.0 0.9 3 0.8 3 0.8 1 1 Fig. 4 | Relaxation maps of the stable TPD glass showing the impact of Data points correspond to the total liquid area calculated as explained in the dynamic facilitation. a, Growth of several liquid regions as a function of time at Supplementary Information. The error bars show an uncertainty proportional to Tg1 + 16 K. The growth time is normalized with respect to the observed nucleation the total length of the perimeter of the liquid regions, and are thus considered to time for each nucleus. The radius of the SCL regions represent best fit values to be 10% of the transformed fraction. c,d, Spatio-temporal maps of the relaxation the AFM images. The error bars range between the smallest and largest circles of the glass including lateral progression of the liquid by dynamic facilitation that could fit at least some of the points of the perimeter of the liquid regions, as at 34% (c) and 90% (d) of transformed fraction. The colour scale marks the time shown in Supplementary Fig. 12. b, Transformed fraction of liquid versus time. evolution. The nuclei evaluated in a are indicated by numbers. the maps at 90% of the transformation of Figs. 3b and 4d evidence highly stable, inactive regions that may only transform from adjacent that most of the glass transforms into liquid by dynamic facilitation, zones that have become liquid. Roughly 10–15% of the glass has still not suggesting that a large fraction of the as-prepared glass consists of transformed after 1.5 × 104 s at Tg1 + 16 K (darker brown areas in Fig. 4d), Nature Physics | Volume 19 | October 2023 | 1509–1515 1513 Article https://doi.org/10.1038/s41567-023-02125-0 indicating that these regions show relaxation times above roughly 8. Cicerone, M. T., Zhong, Q. & Tyagi, M. Picosecond dynamic 2 × 106τα. The mechanism described above is consistent across a range heterogeneity, hopping, and Johari-Goldstein relaxation in of annealing temperatures above Tg1 and thickness of the TPD layers glass-forming liquids. Phys. Rev. Lett. 113, 117801 (2014). (60–400 nm), as briefly shown in Supplementary Figs. 7–9 and 11, but 9. Caporaletti, F. et al. Experimental evidence of mosaic structure in we leave a complete study for future work. We emphasize trilayers strongly supercooled molecular liquids. Nat. Commun. 12, 1867 grown on clean and smooth surfaces transform by a ‘bulk-like’ process (2021). that is not influenced by interfaces in the thickness range analysed here, 10. Paeng, K., Park, H., Hoang, D. T. & Kaufman, L. J. Ideal probe in agreement with previous measurements33,43. The onset temperature single-molecule experiments reveal the intrinsic dynamic of formation of the liquid regions is the same, Ton = 79.5 ± 0.5 °C, for both heterogeneity of a supercooled liquid. Proc. Natl Acad. Sci. USA 13/63/13 and 13/200/13 trilayers, as demonstrated during in operando 112, 4952–4957 (2015). AFM and/or optical microscopy measurements under slow heating rates 11. Zhang, P., Maldonis, J. J., Liu, Z., Schroers, J. & Voyles, P. M. of 0.1 K min−1 (Supplementary Fig. 9). This value also agrees with the Spatially heterogeneous dynamics in a metallic glass forming Ton estimated from differential scanning calorimetry data at 0.1 K min−1 liquid imaged by electron correlation microscopy. Nat. Commun. (Ton = 79.5 ± 0.5 °C) on a bulk, 20-µm-thick sample vapour-deposited at 9, 1129 (2018). 0.85Tg. (Supplementary Fig. 9, right). We also note in passing that the 12. Wisitsorasak, A. & Wolynes, P. G. Dynamical heterogeneity of the dynamics can be accelerated by the presence of defects or contamina- glassy state. J. Phys. Chem. B 118, 7835–7847 (2014). tion on the substrate surface before vapour deposition (or increased 13. Scalliet, C., Guiselin, B. & Berthier, L. Thirty milliseconds in the life roughness of the substrate) as shown in Supplementary Figs. 10 and 11. of a supercooled liquid. Phys. Rev. X 12, 041028 (2022). 14. Vila-Costa, A., Gonzalez-Silveira, M., Rodríguez-Tinoco, C., Conclusions and outlook Rodríguez-López, M. & Rodriguez-Viejo, J. Emergence of We have reported a clear, direct and real-time visualization of the equilibrated liquid regions within the glass. Nat. Phys. 19, 114–119 spatio-temporal dynamics of the devitrification of an ultrastable glass (2023). by AFM and optical measurements. Our data unambiguously dem- 15. Kearns, K. L., Ediger, M. D., Huth, H. & Schick, C. One micrometer onstrate the correlation between time and length scales in glasses. length scale controls kinetic stability of low-energy glasses. Glasses with ultraslow average dynamics and huge transformation J. Phys. Chem. Lett. 1, 388–392 (2010). times, τglass ≅ 106τα, show giant length scales, measured as the aver- 16. Rodríguez-Tinoco, C. et al. Surface-bulk interplay in age distance between patches of fast mobility that equilibrate into vapor-deposited glasses: crossover length and the origin of front the liquid, of the order of 100 μm at Tg1 + 16 K. This giant length scale transformglasses with exceptional thermodynamic and kinetic paves the way to manipulate the properties of ultrastable glasses by stabilityation. Phys. Rev. Lett. 123, 155501 (2019). nanostructuration, opening an interesting avenue to further influence 17. Vila-Costa, A. et al. Nucleation and growth of the supercooled the dynamics. We also anticipate that reducing the stability of the glass liquid phase control glass transition in bulk ultrastable glasses. (that is, its time scale) will also result in a reduction of the length scale. Phys. Rev. Lett. 124, 076002 (2020). By using the mechanical instabilities described here, it should be pos- 18. Swallen, S. F. et al. Organic glasses with exceptional sible to test the existence of heterogeneity on a smaller time and spatial thermodynamic and kinetic stability. Science 315, 353–356 scale by using ultrathin capped layers in less-stable glasses, ultimately (2007). enabling a quasi-direct observation of dynamic heterogeneities in 19. Ediger, M. D. Perspective: highly stable vapor-deposited glasses. liquid-cooled glasses. J. Chem. Phys. 147, 210901 (2017). 20. Rodriguez-Tinoco, C., Gonzalez-Silveira, M., Ramos, M. A. & Online content Rodriguez-Viejo, J. Ultrastable glasses: new perspectives for an Any methods, additional references, Nature Portfolio reporting sum- old problem. La Riv. del. Nuovo Cim. 45, 325–406 (2022). maries, source data, extended data, supplementary information, 21. Queen, D. R., Liu, X., Karel, J., Metcalf, T. H. & Hellman, F. Excess acknowledgements, peer review information; details of author contri- specific heat in evaporated amorphous silicon. Phys. Rev. Lett. butions and competing interests; and statements of data and code avail- 110, 135901 (2013). ability are available at https://doi.org/10.1038/s41567-023-02125-0. 22. Herrero, C., Scalliet, C., Ediger, M. D. & Berthier, L. Two-step devitrification of ultrastable glasses. Proc. Natl Acad. Sci. U.S.A. References 120, e2220824120 (2023). 1. Ediger, M. D. Spatially heterogeneous dynamics in supercooled 23. Bowden, N., Brittain, S., Evans, A. G., Hutchinson, J. W. & liquids. Annu. Rev. Phys. Chem. 51, 99–128 (2000). Whitesides, G. M. Spontaneous formation of ordered structures in 2. Karmakar, S., Dasgupta, C. & Sastry, S. Growing length scales and thin films of metals supported on an elastomeric polymer. Nature their relation to timescales in glass-forming liquids. Annu. Rev. 393, 146–149 (1998). Condens. Matter Phys. 5, 255–284 (2014). 24. Chung, J. Y., Nolte, A. J. & Stafford, C. M. Surface wrinkling: a 3. Sillescu, H. Heterogeneity at the glass transition: a review. J. Non versatile platform for measuring thin-film properties. Adv. Mater. Cryst. Solids 243, 81–108 (1999). 23, 349–368 (2011). 4. Richert, R. Heterogeneous dynamics in liquids: fluctuations in 25. Cerda, E. & Mahadevan, L. Geometry and physics of wrinkling. space and time. J. Phys. Condens. Matter 14, 201 (2002). Phys. Rev. Lett. 90, 074302 (2003). 5. Vidal Russell, E. & Israeloff, N. E. Direct observation of molecular 26. Yin, J., Cao, Z., Li, C., Sheinman, I. & Chen, X. Stress-driven cooperativity near the glass transition. Nature 408, 695–698 buckling patterns in spheroidal core/shell structures. Proc. Natl (2000). Acad. Sci. USA 105, 19132–19135 (2008). 6. Schiener, B., Böhmer, R., Loidl, A. & Chamberlin, R. V. 27. Hoinka, N. M., Ostwald, C. & Fuhrmann-Lieker, T. Two-dimensional Nonresonant spectral hole burning in the slow dielectric wrinkle resonators for random lasing in organic glasses. Sci. Rep. response of supercooled liquids. Science 274, 752–754 (1996). 10, 2434 (2020). 7. Zondervan, R., Kulzer, F., Berkhout, G. C. G. & Orrit, M. Local 28. Torres, J. M., Bakken, N., Li, J. & Vogt, B. D. Substrate temperature viscosity of supercooled glycerol near Tg probed by rotational to control moduli and water uptake in thin films of vapor diffusion of ensembles and single dye molecules. Proc. Natl Acad. deposited N,N′-di(1-naphthyl)-N,N′-diphenyl-(1,1′-biphenyl)- Sci. USA 104, 12628–12633 (2007). 4,4′-diamine (NPD). J. Phys. Chem. B. 119, 11928–11934 (2015). Nature Physics | Volume 19 | October 2023 | 1509–1515 1514 Article https://doi.org/10.1038/s41567-023-02125-0 29. Bakken, N., Torres, J. M., Li, J. & Vogt, B. D. Thickness 39. Das, A., Banerji, A. & Mukherjee, R. Programming feature size dependent modulus of vacuum deposited organic molecular in the thermal wrinkling of metal polymer bilayer by modulating glasses for organic electronics applications. Soft Matter 7, substrate viscoelasticity. ACS Appl. Mater. Interfaces 9, 7269 (2011). 23255–23262 (2017). 30. Moon, M.-W. et al. Wrinkled hard skins on polymers created 40. Dalal, S. S., Walters, D. M., Lyubimov, I., de Pablo, J. J. & Ediger, M. by focused ion beam. Proc. Natl Acad. Sci. USA 104, 1130–1133 D. Tunable molecular orientation and elevated thermal stability (2007). of vapor-deposited organic semiconductors. Proc. Natl Acad. Sci. 31. Chung, J. Y., Nolte, A. J. & Stafford, C. M. Diffusion-controlled, USA 112, 4227–4232 (2015). self-organized growth of symmetric wrinkling patterns. 41. Jin, Y. et al. Glasses denser than the supercooled liquid. Proc. Natl Adv. Mater. 21, 1358–1362 (2009). Acad. Sci. USA 118, 1–7 (2021). 32. Mondal, K., Liu, Y., Shay, T., Genzer, J. & Dickey, M. D. Application 42. Jack, R. L. & Berthier, L. The melting of stable glasses is governed of a laser cutter to pattern wrinkles on polymer films. ACS Appl. by nucleation-and-growth dynamics. J. Chem. Phys. 144, Polym. Mater. 2, 1848–1855 (2020). 244506 (2016). 33. Ràfols-Ribé, J. et al. Kinetic arrest of front transformation to 43. Zhang, Y. et al. Long-range correlated dynamics in ultra-thin gain access to the bulk glass transition in ultrathin films of molecular glass films. J. Chem. Phys. 145, 114502 (2016). vapour-deposited glasses. Phys. Chem. Chem. Phys. 20, 29989–29995 (2018). Publisher’s note Springer Nature remains neutral with regard to 34. Vila-Costa, A., Gonzalez-Silveira, M., Rodríguez-Tinoco, C., jurisdictional claims in published maps and institutional affiliations. Rodríguez-López, M. & Rodriguez-Viejo, J. Emergence of equilibrated liquid regions within the glass. Nat. Phys. 19, Open Access This article is licensed under a Creative Commons 114–119 (2023). Attribution 4.0 International License, which permits use, sharing, 35. Sepúlveda, A., Swallen, S. F. & Ediger, M. D. Manipulating the adaptation, distribution and reproduction in any medium or format, properties of stable organic glasses using kinetic facilitation. as long as you give appropriate credit to the original author(s) and the J. Chem. Phys. 138, 12A517 (2013). source, provide a link to the Creative Commons license, and indicate 36. Mei, B., Zhuang, B., Lu, Y., An, L. & Wang, Z.-G. Local-average free if changes were made. The images or other third party material in this volume correlates with dynamics in glass formers. J. Phys. Chem. article are included in the article’s Creative Commons license, unless Lett. 13, 3957–3964 (2022). indicated otherwise in a credit line to the material. If material is not 37. Gao, S. & Simon, S. L. Measurement of the limiting fictive included in the article’s Creative Commons license and your intended temperature over five decades of cooling and heating rates. use is not permitted by statutory regulation or exceeds the permitted Thermochim. Acta 603, 123–127 (2015). use, you will need to obtain permission directly from the copyright 38. Rodríguez-Tinoco, C., Rams-Baron, M., Rodríguez-Viejo, J. & holder. To view a copy of this license, visit http://creativecommons. Paluch, M. Emergence of a substrate-temperature-dependent org/licenses/by/4.0/. dielectric process in a prototypical vapor deposited hole-transport glass. Sci. Rep. 8, 1380 (2018). © The Author(s) 2023 Nature Physics | Volume 19 | October 2023 | 1509–1515 1515 Article https://doi.org/10.1038/s41567-023-02125-0 Methods temperature control set to room temperature and after stabilization Growth of organic glasses the resonance frequency was tuned, the laser position corrected and Previous experiments have shown that capping a low Tg stable glass a few images at room temperature were acquired. at both interfaces with a higher Tg material inhibits the formation of With the selected 512 points per line and two lines per second, liquid fronts at its surfaces and interfaces forcing the glass to devitrify each image took 256 s to complete (more than 4 min). Lower acquisi- on heating through a bulk-like mechanism16,17,33,35,44. The stack we use tion times were tested resulting in poor resolution images. Although in this work consists of trilayers of TCTA–TPD–TCTA grown on a rigid shorter times would be highly desirable, the evolution of the features substrate (native SiO2–Si) or bilayers of TCTA–TPD–rigid substrate was slow enough to be tracked satisfactorily. During withdrawal and with various TPD and TCTA thicknesses. TCTA has a glass transition engage cycles, the tip would not land exactly on the same position (due temperature of 428 K (Tg2) and the Tg of TPD is 333 K (Tg1), so annealing to hysteresis of the piezoscanner), which makes it difficult to localize above 333 K induces the transformation of the TPD while the TCTA a specific region in transformed samples. remains as a rigid layer. At the target temperature, the apparent lateral dimensions (for Thin-film glasses of TPD and TCTA were grown on both Si sub- example, periodicity of the wrinkles) become smaller due to the strates and SiNx membranes at a substrate temperature of 0.85Tg1 for increase of the strain coefficients (nm V−1) of the piezoscanner with TPD and 0.67Tg2 for TCTA (Tsub = 285 K), to ensure TPD is formed as temperature that is not corrected during image acquisition51. The ultrastable glass. Since the TPD thickness is only around 65 nm, it is actual dimensions should be increased by a factor less than 10% (typi- capped at both sides by an ultrathin layer of a TCTA layer of 13 nm. The cally 5%), as obtained from images acquired at room temperature after growth rate is fixed at 0.08 nm s−1 for both materials. The trilayer struc- the isothermal treatment. A similar trend is expected for the heights ture is schematically shown in Fig. 1a (left). The kinetic and thermody- but in this case the correction is more involved since the tip–surface namic stability of the TPD layers is measured before and after capping interaction should change with temperature, a behaviour that has not with ultrathin layers of TCTA by in situ fast-scanning nanocalorimetry, been explored sufficiently here, combined with the difficulty to localize as already shown elsewhere44–46. a specific region at room temperature that was previously measured at high temperature, so that the obtained heights remain as apparent In operando AFM (uncorrected), although a 10% increase would be reasonable51. Further Temperature-dependent AFM experiments were performed with a calibration experiments should be undertaken to quantify the variation Keysight Technologies 5500 AFM controlled with a PicoScan 3000 of height with temperature. The spatial periodicity of the transformed electronics system, and enclosed in a chamber that was vibrationally surfaces has been calculated from the power spectral density (PSD) and acoustically isolated47. Relative humidity was reduced below 10% distribution as a function of the wavevector modulus49. In this case, relative humidity by circulating nitrogen. Heating was performed by the wavelength corresponds to the inverse of the wavevector (without gluing the samples with silver paste to a modified commercial high the 2 π factor). temperature stage from Keysight Technologies. The temperature was controlled with a LakeShore 311 Temperature Controller in a closed loop Finite element simulations using the internal PT-100 sensor of the heater and measured with an The thin-film system was modelled as an axisymmetric disc stiffened additional PT-100 sensor glued with silver paste to the surface of the by the semi-infinite Si substrate (Supplementary Fig. 3). The overall copper block of the heater (and close to the sample) and the resistance thin-film geometry was discretized using 20-noded, brick-shape ele- measured in a four-wire configuration. The temperature at the sample ments. Since possible processing stresses are unknown, the disc was surface was periodically calibrated with an additional PT-100 sensor assumed to be in the fully relaxed (or stress-free) configuration at the glued with silver paste on a silicon wafer of the same thickness as that deposition temperature of 285 K. Linear thermoelasticity was then of the samples. The x, y and z dimensions of the piezoscanner were invoked as the temperature was ramped to 349 K in the first load step of calibrated with reference samples. For the 1 × 1 μm2 scans, PS-PMMA the finite element simulations. This temperature was set in accordance block copolymers with fingerprint-like morphology and periods of with that where the surface corrugation occurred in the experiments. 35 nm were used48, while for larger scans (up to 20 × 20 μm2) commer- An inner disc of radius r representing the SCL region was then weak- cial calibration grating (TGXYZ01 and TGXYZ02 from MikroMasch) ened in a second simulation stage, where numerical convergence was were used. Images were analysed with both WSxM49 and Gwyddion50. ensured after a slight temperature increment of 0.05 K. The finite ele- The AFM was operated in open-loop intermittent contact mode ment simulations were performed using ANSYS Multiphysics (release using microfabricated silicon cantilevers with force constants of about v.19.0, http://www.ansys.com). 42 N m−1 and ultrasharp silicon tips (nominal tip radius R < 10 nm) The cylindrical, weakened SCL regions were assumed to deform in (PPP-NCHR, Nanosensors). The resonance frequencies were in the accordance with linear elasticity with specific Young’s modulus E and 250–350 kHz range. After laser alignment, the sample was fixed to the Poisson’s ratio ν. Alternatively, the neo-Hookean hyperelastic formalism heater and the heater stage coupled to the AFM head. When the whole was also invoked so as to account for the potential decrease of the elastic setup was vibrationally and acoustically isolated, several images rang- modulus of the weakened SCL regions at small strains. In this model, 2 ing from 1 × 1 up to 20 × 20 μm2 in scan size were taken until stabilization the total strain energy density is given by U = C(I1 − 3) + d−1 (detF − 1) , of the microscope. Once stabilized at room temperature, the sample where C and d are material parameters, I1 is the first deviatoric strain was withdrawn and the heater was driven up to Tg1 −10 K at a heating rate invariant and detF is the volume ratio. The second term in the right-hand of 5 K min−1 (as recommended by the manufacturer). Once the desired side of the neo-Hookean formulation explicitly accounts for possible temperature was reached, the resonance frequency was tuned (the material compressibility or volume changes, so that the SCL regions frequency decreases for increasing temperature), the laser was aligned approach full incompressibility when d → 0. again, the tip was engaged and a few images were acquired to check The capping TCTA top and bottom layers were assumed to remain for the stability of the instrument. Then, the sample was withdrawn linear elastic in all simulations, with matching Young’s modulus E and again and the temperature increased up to the target temperature, Poisson ratio ν to those in the glassy state (see the thermoelastic prop- Tg1 + 16 K. Again, when the temperature was reached, the tuning of the erties for the thin-film materials and Si substrate in Supplementary resonance frequency was not necessary in most cases to adjust the laser Information). The material parameters for the SCL disc where then position. Then, the tip was engaged to find suitable regions to follow allowed to vary within the wider ranges of C = 3.84 × 106 to 1.43 × 107 Pa the isothermal evolution of the surface as a function of time. When the and d = 1.12 × 10−7 to 2.79 × 10−11Pa−1. The experimentally observed transi- temperature experiments were finished, the sample was withdrawn, the tion between the onset of a single ridge to wrinkles with increasing r Nature Physics Article https://doi.org/10.1038/s41567-023-02125-0 from 250 to 1,000 nm was then mimicked within the simulation subset Acknowledgements where the SCL properties varied as before. Grant nos. PID2020-117409RB-I00 (J.R.-V., C.R.-T. and M.G.-S.), PGC2018-095032-B-100 (J.F.), PID2020-115663GB-C3 Data availability (J.A.P.) and PID2019-106744GB-I00 (J.A.), funded by MCIN/ Data that support the plots and other findings of this study are available AEI/10.13039/501100011033. Project 2021SGR-00644 funded at https://doi.org/10.5281/zenodo.7950765. by AGAUR. C.R.-T. is a Serra Hunter Fellow. The ICN2 is funded by the CERCA programme/Generalitat de Catalunya. The ICN2 is References supported by the Severo Ochoa Centres of Excellence programme, 44. Leon-Gutierrez, E., Sepúlveda, A., Garcia, G., Clavaguera-Mora, grant no. CEX2021-001214-S, funded by grant no. MCIN/ M. T. & Rodríguez-Viejo, J. Stability of thin film glasses of toluene AEI/10.13039.501100011033. and ethylbenzene formed by vapor deposition: an in situ nanocalorimetric study. Phys. Chem. Chem. Phys. 12, Author contributions 14693–14698 (2010). J.R.-V., C.R.-T. and M.G.-S. conceived the project. J.F., M.R.-R. and 45. Leon-Gutierrez, E., Sepúlveda, A., Garcia, G., Clavaguera-Mora, M. A.V.-C. performed the AFM experiments. T.B. performed the optical T. & Rodríguez-Viejo, J. Correction: stability of thin film glasses of microscopy measurements. T.B., J.F. and M.R.-R. analysed the AFM toluene and ethylbenzene formed by vapor deposition: an in situ and optical microscopy images. J.A.P. and J.A. contributed with the nanocalorimetric study. Phys. Chem. Chem. Phys. 18, 8244–8245 mechanical analysis and interpretation. C.R.-T., M.G.-S. and J.R.-V. (2016). discussed fundamental ideas of the transformation mechanisms. 46. Rodríguez-Tinoco, C. et al. Evaluation of growth front velocity J.R.-V. wrote the manuscript with C.R.-T., J.F., J.A., J.A.P. and M.G.-S. All in ultrastable glasses of indomethacin over a wide temperature authors commented on the manuscript. J.R.-V. supervised the project. interval. J. Phys. Chem. B 118, 10795–10801 (2014). 47. Verdaguer, A., Cardellach, M. & Fraxedas, J. Thin water films Competing interests grown at ambient conditions on BaF2(111) studied by scanning The authors declare no competing interests. polarization force microscopy. J. Chem. Phys. 129, 174705 (2008). Additional information 48. Fernández-Regúlez, M., Evangelio, L., Lorenzoni, M., Fraxedas, Supplementary information The online version contains supplementary J. & Pérez-Murano, F. Sub-10 nm resistless nanolithography for material available at https://doi.org/10.1038/s41567-023-02125-0. directed self-assembly of block copolymers. ACS Appl. Mater. Interfaces 6, 21596–21602 (2014). Correspondence and requests for materials should be addressed to 49. Horcas, I. et al. A software for scanning probe microscopy Cristian Rodríguez-Tinoco or Javier Rodriguez-Viejo. and a tool for nanotechnology. Rev. Sci. Instrum. 78, 013705 (2007). Peer review information Nature Physics thanks Heinz Sturm and the 50. Nečas, D. & Klapetek, P. Gwyddion: an open-source software for other, anonymous, reviewer(s) for their contribution to the peer review SPM data analysis. Open Phys. 10, 181–188 (2012). of this work. 51. Ma, C., Chen, Y. & Huang, W. Environmental temperature effect on dimensional measurements of atomic force microscopy. Reprints and permissions information is available at Nanotechnol. Precis. Eng. 4, 023003 (2021). www.nature.com/reprints. Nature Physics