Ultrafast electron crystallography reveals the atomic pathway of a light-driven correlated insulator-to-metal transition

Ultrashort laser pulses can transiently reshape the electronic, magnetic, and structural properties of materials, in some cases creating states that cannot be stabilized under equilibrium conditions [1, 2, 3]. Such capabilities have opened routes to ideas such as light-induced superconductivity [4], all-optical switching of magnetic order on picosecond timescales [5], and metastable structural phases that underpin commercial phase-change memory [6]. Nonetheless, determining the atomic structures and transition pathways of these short-lived states has remained beyond the reach of full crystallographic refinement from ultrafast diffraction data. The limitation is intrinsic: obtaining complete, high-precision diffraction datasets on ultrafast timescales is challenging because transient scattering signals are weak and measurements must balance temporal resolution, reciprocal-space coverage, acquisition time, and sample damage. As a result, time-resolved electron and x-ray diffraction experiments often track one or a few selected selected Bragg reflection intensities or structural order parameters [7, 8, 9]. This contrasts with static electron and x-ray diffraction captured under equilibrium conditions, where extensive datasets of three-dimensional diffraction have enabled precise determination of fractional atomic coordinates within the unit cell [10, 11].

In this work, we overcome this limitation by extending three-dimensional electron diffraction (3DED) into the ultrafast regime, a technique we call ultrafast 3DED [12]. Static 3DED is an established method for atomic-scale structure determination of sub-micron crystals across materials science and structural biology, where single-crystal x-ray diffraction is often not feasible [13, 14]. Its strength lies in measuring Bragg intensities over a series of crystal tilt angles, thereby sampling the crystal from multiple orientations relative to the electron beam. We combine this tilt-series acquisition with an optical pump-electron probe scheme (Figure 1a), accumulating diffraction signal stroboscopically over many pump-probe cycles and merging Bragg intensities measured over the full tilt series at each pump-probe delay. The resulting 3DED dataset is refined directly to obtain the transient atomic coordinates within the unit cell during the photoexcited phase transition (Figure 1b). By acquiring diffraction data off the main-zone axes and integrating Bragg reflections over multiple tilt angles, ultrafast 3DED suppresses multiple scattering and improves intensity accuracy relative to conventional 2D diffraction measurements, allowing the average transient crystal structure to be determined at every pump–probe delay rather than inferred from changes in selected Bragg reflections.

We demonstrate 3DED on vanadium dioxide ($\mathrm{VO}_{2}$), a prototypical strongly correlated material, in which the atomic pathway of the photoinduced insulator-to-metal transition has remained unresolved despite decades of study. This material undergoes a coupled structural and electronic phase transition near 340 K, from a monoclinic (M1) insulator with dimerized vanadium pairs to a tetragonal metallic rutile (R) phase (Figure 1c). The M1 distortion relative to the rutile parent can be decomposed into three coupled structural components: (i) a monoclinic shear, reflected in the deviation of the monoclinic angle $\beta$ (defined between the $a$ and $c$ lattice vectors) from the rutile-reference value of approximately 116° in the same monoclinic setting, corresponding to what is sometimes described as unit-cell distortion; (ii) V–V dimerization along the chains of vanadium atoms running parallel to the rutile $c_{\mathrm{R}}$ axis, in which neighboring V atoms form alternating short intra-dimer and long interdimer V–V distances; (iii) a lateral zigzag displacement, in which the V atoms shift transversely to the chain in an alternating pattern, corresponding to what is sometimes described as dimer tilting or chain twisting. These components vanish in the ideal rutile structure and together define the monoclinic M1 distortion [7, 8, 15]. This decomposition provides a structural framework for describing the photoinduced M1-to-R pathway and for identifying the components of the M1 distortion driving the electronic gap opening.

Recent experimental and theoretical studies suggest that the M1-to-R transition may proceed through pathways more complex than a direct structural interpolation. Here, a direct, linear pathway denotes a coordinated evolution in which the atomic coordinates move along a straight line in three-dimensional space between the M1 and R structures. A nonlinear or asynchronous pathway instead denotes a trajectory that deviates from this straight structural interpolation, either because different components evolve on different timescales or because the system passes through an intermediate crystallographic structure. Li et al. [8] reported that V–V dimerization and lateral zigzag initially evolve toward their rutile values along a linear, common-timescale trajectory. At later delays, however, the motion no longer follows this trajectory: the lateral zigzag continues to relax, whereas the V–V dimerization changes comparatively little. Baum et al. [7] inferred a stepwise, temporally separated pathway: a femtosecond component assigned to V–V bond dilation, followed by picosecond atomic displacements and slower shear motion. Johnson et al. [9] found that incoherent, spatially correlated lattice distortions, identified as interacting polarons, create a transient state in which the energy barrier to the metallic phase is reduced. Theoretical work by Grandi et al. [16] proposed a two-step transition in which V–V dimerization is lost before the lateral zigzag distortion, producing an intermediate monoclinic metallic phase between the M1 insulator and the rutile metal.

While these works significantly advanced our understanding of the mechanisms behind the ultrafast phase transition, they do not determine the complete transient crystallographic structure; instead, the transition pathway is inferred indirectly from projected coordinates, Bragg-intensity classes, diffuse-scattering correlations, or theoretical order-parameter landscapes. In contrast, ultrafast 3DED refines a structural model against Bragg intensities measured over a series of crystal tilt angles. Each pump–probe delay yields a three-dimensional dataset containing several hundred Bragg reflections, sufficient to refine the fractional coordinates of all symmetry-inequivalent atoms in the unit cell. This capability extends ultrafast electron diffraction beyond order-parameter tracking to time-resolved crystallographic refinement, enabling direct tracking of the atomic trajectories throughout the photoexcited phase transition.

In this work, we excite a structural phase transition in $\mathrm{VO}_2$ using laser pulses with a photon energy of 1.2 eV and a duration of 300 fs. The resulting changes in diffraction are probed using 1.5 ps electron pulses in a stroboscopic approach at a repetition rate of 12 kHz. The sample is a focused ion beam (FIB) lamella prepared from a single-crystalline $\mathrm{VO}_2$ micro-needle [17], with dimensions of $4 \times5\ \mu\mathrm{m}^{2}$ and an approximate thickness of 100 nm (see SM). The experimental data are processed using the REDp software [18], which reconstructs the time-resolved 3D reciprocal lattice. The intensities of the Bragg reflections in each 3DED dataset at a given time delay are integrated using XDS [19]. These diffraction patterns can be projected along arbitrary zone axes (see fig:2a).

To monitor the phase transition, we track the total intensity of the Bragg peaks (see fig:2b), splitting them into two categories: M1-exclusive peaks, which are present in the monoclinic M1 phase but absent in the rutile R phase, and shared peaks, common to both structures. The M1-exclusive peaks are sensitive to component (ii) and (iii) of the M1 distortion. At each delay, the signal is averaged over more than 500 such reflections. The M1-exclusive intensity follows a double-exponential decay (see fig:2b), with both the bi-exponential dynamics and the associated timescales consistent with previous reports [7, 8, 20]: a fast decay reduces the intensity to approximately 50% within $\sim10$ ps, and a slower decay completes the suppression on a characteristic timescale of 400 ps. The intensity loss at the M1-exclusive peaks reflects the decrease of the dimerization amplitude. In the M1 phase, the V atoms are organized into two interleaved sub-lattices displaced from each other by the dimerization, with twice the rutile periodicity along the chain. This doubling of the periodicity along the rutile $c_{\mathrm{R}}$ direction produces the M1-exclusive peaks and suppresses the structure factor at the shared peaks. As the photoexcited transition proceeds and the dimerization decreases, the two sub-lattices merge into a single uniform lattice: the M1-exclusive peaks are suppressed, and the intensity of the shared peaks grows. This anticorrelated trend is observed in our data, consistent with the redistribution of Bragg intensity as the M1 superstructure is suppressed.

Beyond intensity tracking, ultrafast 3DED allows us to determine the unit-cell parameters and the atomic fractional coordinates (the positions of atoms expressed as fractions of the unit-cell dimensions) at each pump-probe delay. We use this capability to track the lattice parameters of the unit cell, as well as the individual atomic coordinates within it. Of particular interest for $\mathrm{VO}_2$ is the time-dependent evolution of the monoclinic angle $\beta$, the angle between the $a$ and $c$ lattice vectors of the unit cell (Fig. 2c), which directly tracks component (i), the monoclinic shear. In contrast to the V–V dimerization and zigzag displacement (Fig. 2d,e), which exhibit a pronounced fast response immediately following photoexcitation, the evolution of $\beta$ is dominated by a slower exponential relaxation. The dimerization and zigzag components also display a slower evolution on a comparable timescale. This indicates that the slow relaxation of components (i)–(iii) toward their rutile values proceeds on a common timescale rather than asynchronously, even though components (ii) and (iii) additionally exhibit a distinct ultrafast response immediately following photoexcitation.

Although our measurements may not accurately capture small absolute changes in $\beta$ because of the limited number of tilt angles, they determine its temporal evolution with high precision and reveal a clear trend. The high precision arises because the experimental conditions, such as tilt angles and laser illumination, remain consistent across all delays. Static measurements report $\beta= 115.2^\circ$ for the monoclinic phase and approximately $116^\circ$ for the rutile phase in the same monoclinic representation [15, 21].

We employ two complementary approaches to determine the transient atomic coordinates during the phase transition. First, we use the program SHELXL [22, 23] to refine the atomic coordinates of the vanadium (Fig. 3a) and two oxygen atoms (see SM) using the full reflection set, with all monoclinic structures successfully converging. This provides the absolute crystallographic coordinates within the refined transient unit cell. Second, using only the M1-exclusive reflections, we apply a least-squares refinement to determine the relative, time-dependent atomic coordinates. This provides an independent determination of the transient M1 structural evolution that excludes reflections that may contain contributions from the rutile phase. Both approaches yield excellent agreement of the V coordinates (3), providing a consistency check on the extracted transient structural evolution.

The refined V fractional coordinates capture both components (ii) and (iii) of the M1 distortion: the longitudinal V–V pairing and the transverse zigzag displacement, respectively. As shown in fig:3a, both structural components evolve from their M1 values toward the rutile configuration with a common overall relaxation behavior. The most abrupt changes occur within the first few picoseconds, during which the atoms traverse approximately half of the total structural displacement. This initial motion is followed by a much slower evolution extending beyond 1 ns, during which the remaining displacement is gradually completed. The vanadium and oxygen coordinates therefore follow nearly straight trajectories from their M1 to their R positions, defining a direct, linear pathway in coordinate space, as illustrated by the correlation plot in fig:3b. Additional correlation plots of the vanadium fractional coordinates are shown in fig:S7 in the SM. The $x$-components of the fractional coordinates of the O1 and O2 atoms are shown in fig:3c,fig:3d, respectively. All remaining components and correlation plots are shown in fig:S2. A residual V–V dimerization of approximately 15% of its initial M1-to-R distortion nevertheless persists into the nanosecond regime, indicating that a small fraction of the initial V–V dimerization displacement remains despite the average crystallographic structure being indexed as rutile. Potential phase heterogeneity [24] is not expected to significantly affect the determination of the M1 atomic positions. Fractional coordinates obtained from the least-squares refinement to M1-exclusive reflections (see SM) are in close agreement with those obtained using the full reflection set in SHELXL, which also includes reflections potentially shared with the rutile phase. This agreement indicates that inclusion of the potentially phase-overlapping reflections does not significantly affect the refined transient M1 fractional coordinates.

Having tracked all three components of the M1 distortion as a function of pump–probe delay, we next identify which structural distortion is responsible for the collapse of the electronic band gap. The insulating M1 phase arises from the interplay of structural symmetry breaking and electronic correlations. In particular, the formation of V–V dimers along the rutile chains splits the chain-aligned $d_{\parallel}$ orbitals into bonding and antibonding states, while the accompanying distortion of the VO$_6$ octahedra shifts the relative energy of the $e_g^{\pi}$ states (fig:4a). Together, these distortions create the structural framework for the insulating M1 phase, while electronic correlations determine the magnitude and stability of the gap. The relative importance of the individual structural distortions, however, has remained under debate.

The three structural components affect the electronic states of VO$_2$ in different ways. The monoclinic shear changes the unit-cell geometry and modifies the V–O–V bond angles, thereby altering the dispersion of the electronic bands. However, this distortion alone does not directly split the V-derived states responsible for the insulating character of the M1 phase. In contrast, V–V dimerization along the rutile chains creates alternating short and long V–V bonds. This changes the overlap between neighboring V orbitals and splits the chain-aligned $d_{\parallel}$ states into bonding and antibonding states, providing the structural basis for the electronic gap formation. The lateral zigzag distortion primarily modifies V–O hybridization and shifts the relative energy of the $e_{g}^{\pi}$ states, thereby influencing the size of the gap, but does not generate the bonding–antibonding splitting of the $d_{\parallel}$ manifold. These considerations suggest that V–V dimerization is the structural distortion most directly linked to the electronic gap opening, consistent with the Goodenough–Eyert picture [25, 26, 27]. Nevertheless, determining the quantitative contribution of each distortion requires explicit calculations of the correlated electronic structure, as the M1 gap emerges from the interplay of structural changes and electronic interactions.

To quantify how the experimentally observed structural changes affect the electronic properties, we perform cluster dynamical mean-field theory (DMFT) calculations on the experimentally refined transient structures. The correlated subspace is constructed from V-centered $t_{2g}$ orbitals, with each V–V dimer treated as a correlated impurity cluster. The calculations include local Coulomb interactions and Hund coupling, together with a nonlocal exchange contribution between the two V atoms of the dimer (see SM). The impurity problem is solved using continuous-time hybridization-expansion quantum Monte Carlo (CT-HYB) [28, 29, 30] within the all-electron DFT+DMFT framework implemented in Relativistic Spin Polarized toolkit (RSPt) [31]. Calculations are performed at $T = 300$ K, corresponding to the experimental conditions. For the equilibrium M1 structure, this approach yields an insulating gap of 0.67 eV, in excellent agreement with the experimentally reported gap of approximately 0.6–0.7 eV from optical spectroscopy, photoemission, and transport measurements [32, 33, 34, 35] (Figure 4b). The cluster-DMFT solution is paramagnetic, and does not impose static long-range magnetic order, while capturing singlet correlations within the V–V dimers. This is consistent with the absence of long-range magnetic order in the M1 phase at room temperature [36, 37].

We next use the experimentally refined M1 structure as a reference to disentangle the structural and many-body contributions to the $d_{\parallel}$ bonding–antibonding splitting. The calculated splitting is approximately 3 eV (Figure 4b). V–V dimerization provides the one-body contribution by creating contrasting hopping amplitudes across the short and long V–V bonds, while many-body effects are described by the full frequency-dependent intersite self-energy of the dimer. In cluster DMFT, this self-energy arises from the on-site interactions $U$ and $J$ together with an explicit intersite exchange interaction $V_{\mathrm{ex}}$, which is included as a full two-body term in the cluster impurity interaction tensor (see SM). At a fixed frequency, the intersite hopping and intersite self-energy enter with opposite signs for the bonding and antibonding combinations, thereby modifying their separation, whereas site-local terms contribute equally to both. With $V_{\mathrm{ex}}=-0.49$ eV together with $U$ and $J$, the calculated gap is 0.67 eV and the occupied $d_{\parallel}$ bonding peak lies 0.6 eV below the chemical potential. Retaining $U$ and $J$ but setting $V_{\mathrm{ex}}=0$, with all other parameters unchanged, the system remains insulating, but the gap decreases to approximately 0.14 eV and the bonding peak shifts to only 0.2 eV below the chemical potential. We then completely remove the interaction tensor, setting $U=J=V_{\mathrm{ex}}=0$, while retaining the same refined M1 geometry. In this one-body limit, the $d_{\parallel}$ splitting decreases to approximately 1.5 eV and the separation between the occupied $d_{\parallel}$ bonding states and the low-lying $e_g^{\pi}$ states is lost (Figure 4c). Together with the shear-only control below, these results identify V–V dimerization as the dominant structural origin of the bonding–antibonding splitting, while showing that many-body interactions substantially enhance this splitting and stabilize the full insulating M1 state.

At the Hartree–Fock level, the first-order off-diagonal self-energy generated by $V_{\mathrm{ex}}$ maps onto the static intersite correction of DFT+$V$, giving $V=-2V_{\mathrm{ex}}=0.98$ eV; the factor of two is the spin sum in the exchange contraction [38] (see SM). Cluster DMFT is not restricted to this static limit: $V_{\mathrm{ex}}$ is sampled as a full two-body interaction, including opposite-spin spin-exchange terms, and therefore also contributes to the frequency-dependent many-body self-energy. The short V–V bond length, in turn, controls the orbital overlap between neighboring V atoms and therefore strongly influences both the one-body hopping and the many-body intersite self-energy. The hopping across the short bond decreases as the bond lengthens. Microscopically, the screened exchange interaction represented by $V_{\mathrm{ex}}$ is likewise expected to weaken with increasing V–V distance. In a simple atomic-orbital picture, the overlap decreases approximately as $e^{-r/L}$, where $r$ is the short V–V bond length and $L$ is an effective orbital decay length. As the exchange matrix element contains an orbital-overlap factor at each electron coordinate, $\lvert V_{\mathrm{ex}}\rvert$ is expected to decrease approximately as $e^{-2r/L}$. We verified numerically that the magnitude of the intersite self-energy decreases as $\lvert V_{\mathrm{ex}}\rvert$ is reduced. The structure itself further modifies the intersite self-energy through changes in hopping, hybridization, and the intersite density matrix, which decreases as the short and long V–V bonds equalize. Thus, as V–V dimerization weakens along the measured photoinduced structural trajectory (fig:2d), the lengthening of the short V–V bond reduces the one-body hopping and weakens the many-body intersite contribution to the bonding–antibonding splitting.

We further examine whether the monoclinic shear component alone can induce an insulating state within the noninteracting one-body limit, $U = J = V_{\mathrm{ex}} = 0$. We construct a hypothetical shear-only structure by placing the rutile fractional coordinates in the monoclinic unit cell, without V–V dimerization or the lateral zigzag displacement. The $d_{\parallel}$ bonding–antibonding splitting is strongly reduced and the calculated spectral function is metallic (fig:4d). Although the shear lifts the symmetry-enforced degeneracy of the $d_{\parallel}$ bands at the rutile Z point (see SM), the resulting splitting is too small to open a gap. This shows that the symmetry breaking associated with the monoclinic cell is not sufficient to produce the insulating M1 electronic structure. For comparison, we also calculate the rutile structure corresponding to the long-delay state, in which all three M1 distortion components vanish. The residual band splitting observed for the shear-only structure then disappears and the $d_{\parallel}$ states recover their rutile symmetry, yielding a metallic spectral function consistent with the experimentally observed metallic steady state at delays beyond 1 ns (fig:2d).

Beyond the interactions within the V–V dimer, V–O hybridization contributes to the insulating gap through both one-body band alignment and ligand-mediated many-body effects. The choice of correlated orbitals determines how these ligand-mediated interactions are incorporated into the many-body treatment. For predominantly V-centered orbitals, an on-site Hubbard interaction alone is insufficient to reproduce the M1 gap: single-site DMFT and spin-polarized DFT+$U$ remain metallic for $U$ values up to 4 eV when the explicit intersite exchange interaction is omitted. When the correlated orbitals contain more O-ligand character, the nominal on-site interaction acts on hybridized V–O orbitals and therefore incorporates contributions that are nonlocal in an atom-centered representation, including screened V–O Coulomb interactions and ligand-mediated contributions to the effective V–V interactions. In this ligand-rich subspace, the $d_{\parallel}$ bonding–antibonding splitting increases and a gap opens even without the explicit intersite exchange interaction. V–O hybridization therefore affects the gap in two complementary ways: ligand-mediated interactions enhance the $d_{\parallel}$ bonding–antibonding splitting, and one-body V–O hybridization controls the relative alignment of the $d_{\parallel}$ and $e_g^{\pi}$ states. The ligand-rich construction, however, does not reproduce a systematic collapse of the gap when dimerization is removed (see SM); the main calculations therefore use the predominantly V-centered subspace.

These results show that the M1 gap cannot be explained by the structural $d_{\parallel}$ bonding–antibonding splitting alone. V–V dimerization provides the essential one-body framework. The full frequency-dependent intersite self-energy, generated by the on-site interactions and the explicit intersite exchange interaction, further increases the $d_{\parallel}$ bonding–antibonding splitting. Ligand-mediated V–O interactions also contribute to this splitting, while one-body V–O hybridization controls the relative alignment of the $d_{\parallel}$ and $e_{g}^{\pi}$ states. Together, these electronic effects stabilize the full insulating state. A complete description of these coupled structural and electronic effects would require a unified treatment of dynamical correlations and nonlocal screening and exchange, such as $\mathit{GW}+\mathrm{DMFT}$ [39, 40, 41, 42]. Recent quasiparticle $\mathit{GW}$ studies of VO$_2$ likewise show that the M1 gap is strongly sensitive to lattice distortions [43].

In conclusion, ultrafast 3DED enables transient crystallography of VO$_2$ by resolving the time-dependent evolution of the unit cell and all symmetry-inequivalent V and O positions. Within the frozen-snapshot equilibrium mapping used here, combining these measurements with cluster DMFT calculations on the experimentally refined transient structures shows that the calculated gap collapse follows primarily from the loss of V–V dimerization along the measured structural trajectory. As the dimerization weakens, the reduced $d_{\parallel}$ overlap diminishes the one-body bonding–antibonding splitting and weakens the structure-dependent intersite self-energy contributions that stabilize the equilibrium M1 phase. Monoclinic shear alone is insufficient to open a gap. V–O hybridization affects both the $d_{\parallel}$ splitting through ligand-mediated interactions and the relative alignment of the $d_{\parallel}$ and $e_{g}^{\pi}$ states. Together, experiment and theory identify M1 VO$_2$ as a correlated dimer insulator and demonstrate how ultrafast 3DED can connect atomic-scale structural dynamics with electronic functionality in correlated materials. The reported trajectory describes the average structure; local distortions invisible to Bragg diffraction may persist.

Acknowledgments

We thank Prof. Bertina Fisher for providing the VO$_2$ microneedles used to create the FIB lamellae in this study. Language polishing was assisted by ChatGPT (OpenAI). All content was critically reviewed and approved by the authors.

Funding: This research was funded by the Knut and Alice Wallenberg Foundation (2012.0321, 2018.0104, 2019.0124), the Swedish Research Council (VR, 2018.00815), and through the ARTEMI national infrastructure (VR 2021-00171, VR 2022-03596 and Strategic Research Council (SSF) RIF21-0026). A.N. acknowledges funding from the Swiss National Science Foundation (SNSF) through project P500PT 214469. H.X. acknowledges funding from the Australian Research Council for the Future Fellowship (FT240100378). OG and CSO acknowledge funding from VR (2025-06120 and 2025-06380). OG, CSO and PT acknowledge funding from European Research Council (ERC), Synergy Grant 854843-FASTCORR. Computational resources were provided by the National Academic Infrastructure for Supercomputing in Sweden (NAISS), funded by the Swedish Research Council. Jianyu Wu gratefully acknowledges the Chinese Scholarship Foundation (CSC) for a doctoral fellowship. L.K. is grateful for the support of the Israel Science Foundation (ISF Grant No. 1397/24).

Competing interests: There are no competing interests to declare.

Data and materials availability: All data are available in the manuscript or the supplementary materials.

References and Notes