Metal-like monoclinic phase and terahertz characteristics in ultrafast phase transition of photoexcited VO₂

Vanadium dioxide (VO $2$⁠) is a typical strongly correlated electron material, which undergoes a reversible phase transition at $T c∼$ 341 K.¹ In the traditional temperature induced phase transition, the relationship of electronic insulator-metal transition (IMT) and structural phase transition (SPT) is coupling, bringing huge controversies on the phase transition mechanism.^2,3 In recent years, photoexcitation has emerged as a powerful way to induce the phase transition of strongly correlated materials.^4–9 The separation of the IMT and the SPT has been theoretically and experimentally discovered in VO $2$ excited by laser pump.^6,8–10 However, the excitation threshold studies^5,11,12 for decoupling the IMT and the SPT, which enables us to reveal the physical mechanism of photoexcited VO $2$ and guide the dynamic control of optoelectronic devices, are still unsolved issues.

In the transition from the monoclinic (M $1$⁠) phase to the rutile (R) phase of VO $2$⁠, its conductivity and the resulting infrared and terahertz (THz) spectral responses dramatically change.^2,3,13,14 Thanks to its notable IMT features, VO $2$ has been widely used in smart sensors,^15,16 rewritable optical memories,^17,18 and dynamically tunable THz devices.^19–22 Especially integrated with functionalized metamaterials,^23,24 various applications in THz regime have been demonstrated. Although the development of femtosecond lasers and terahertz measurement equipment has brought new insights to many experiments, theoretical research on the physical mechanism and terahertz characteristics of photoexcited VO $2$ is still slow. On the one hand, owing to the frequency specificity of terahertz, it is unfeasible to use traditional density functional theory (DFT)^25,26 to calculate terahertz characteristics of crystals by ab initio.²⁷ On the other hand, although the metal-like monoclinic phase structure in the VO $2$ phase transition has been found experimentally,^5,17,28,29 the atomic structure and electronic properties used to achieve THz responses in the ultrafast photoinduced phase transition (PIPT) are still unknown. These reasons hinder the understanding of the physical mechanism in the PIPT, as well as the design and theoretical verification of photoexcited terahertz devices.

In this work, we have investigated the dynamic process of photoexcited VO $2$ by utilizing the Hubbard corrected real-time time-dependent density functional theory (rt-TDDFT $+U$⁠) method^30–32 and discovered the separation with different inducing thresholds for IMT and SPT processes. Under laser irradiations with different intensities, the IMTs occur quasi-instantaneously at 15–45 fs, while the ultrafast SPT processes are completed at 56–141 fs. It can be observed that the metal-like monoclinic phase structures exhibiting the strongest metallicity and the weak insulating monoclinic phases, when the bond lengths of the V–V pairs are closest. In our previous work,²⁷ we successfully calculated the THz characteristics of temperature controlled VO $2$ by using electron–phonon interactions (EPIs) and the equivalent dielectric theory (EMT). This method is extended into the PIPT process in this work. With the EPIs method^33–35 and the Boltzmann transport equations solver,³⁶ the electronic transport properties of metal-like structures are calculated, and, thus, the THz response of the whole photoexcited phase transition can be characterized through the extended quantum-electromagnetic computational method. For verification, the THz transmission spectra of VO $2$ film based the mica substrate are simulated under laser illumination. Great agreement between simulation and measurement validates the correctness of the proposed analysis scheme. The micro and macro characteristics of the whole PIPT demonstrated by our work provide a powerful exploration path and insights to other PIPT materials and the corresponding THz devices.

At the temperature below the transition threshold $∼$341 K, VO $2$ is stabilized in a low-symmetry $M 1$ phase with space group P $2 1/c$ (⁠ $#14$⁠), as presented in Fig. 1(a). As the temperature increases, VO $2$ undergoes first-order phase transition and then transforms into a symmetry R phase with space group P $4 2/mnm$ (⁠ $#136$⁠) shown in Fig. 1(b). The dimerization of V atoms features zigzag V–V chains at the insulating $M 1$ phase, with the theoretical long and short bond lengths $d L=3.153$ Å and $d S=2.489$ Å and the twisting angle $θ= 166.58 °$⁠. At the metallic R phase, the V–V dimerization is broken, resulting in $d R= d L= d S=2.802$ Å and $θ= 180 °$⁠.

Fatbands [Figs. 1(c) and 1(d)] show the contribution of $O 2 p$ and $V 3 d$ orbitals to $M 1$ and R phases, emphasizing that the splitting of $a 1 g$ orbital (⁠ $d x 2 − y 2$⁠) is strongly correlated with the phase transition of VO $2$⁠. The highly directional $a 1 g$ orbital splits into a bonding state $d ∥$ and an antibonding state $d ∥ ∗$⁠, and the shift of $d ∥$ below the Fermi energy leads to a 0.6–0.7 eV bandgap.^37,38 In virtue of this strong electronic localization, the Hubbard correction (⁠ $U$⁠) is involved for the correct calculation of bandgap. The chosen $U=3.1$ eV leads to a 0.68 eV bandgap, in agreement with the theoretical and measured values.^14,27,38 Note that, the $U$ value is selected by numerical search to reduce computational complexity compared with the fist-principles calculation method.²⁷ As the phase transition occurs, $d ∥$ shifts to cross the Fermi energy and couples with $d ∥ ∗$⁠, resulting in the metallic R phase. Therefore, the shift of $d ∥$ is one of the critical issues to investigate the PIPT procedure.

To simulate the photoinduced ultrafast $M 1$-to-R phase transition by the rt-TDDFT + $U$ method, an 800 nm laser pump illuminates the $2×2×2$ supercell of the $M 1$ phase at room temperature (300 K) and non-adiabatic molecular dynamics simulation is performed in real time. In the calculations, the phonon and atomic motion at each time step have been fully considered.

In photoexcited procedure, an external electrical field with the Gaussian waveform, i.e., $E(t)= E 0cos(ωt) e − ( t − t 0 ) 2 / ( 2 σ 2 )$⁠, is employed. Here, the central time $t 0=17.5$ fs and the pulse width $2 σ 2$ = 49 fs. The electric field $E(t)$ is along the a direction of the $M 1$ phase supercell and the electric field strength $E 0$ increases from 0.15 to 0.5 Hartree/bohr. The specific laser signals and the resultant evolution of the excited electron with time are plotted in Fig. S1 in the supplementary material. The lasers excite the valence electrons at the $d ∥$ state to the conduction bands, with the equivalent number of electrons and holes doping. The photoexcited electron (hole) percentages associated with the different electric field amplitudes $E 0$ of the entire laser pulses are presented in Fig. 1(e) and Table S1 in the supplementary material, representing the different laser pump intensities in this simulation work. The hole doping weakens the V–V dimerization and further triggers the IMT and SPT processes.³⁹ 

After carrier doping, the structure dynamic of the $M 1$-to-R transition can be characterized by the mean bond lengths and twisting angles. Figure 2(a) depicts the evolutions of the mean long and short bond lengths (⁠ $d L$ and $d S$⁠) with time below the SPT excitation threshold, wherein the time of $d L$ =  $d S$ is considered as the SPT moment.^10,30 Here, the use of the mean evolutions neutralize the disorder effects from thermal fluctuations. As shown in Fig. 2(a), $d L$ and $d S$ are getting closer as the laser intensity increases from $0%$ to $0.83%$ excited electron. For $0.94%$ excited electron, the complete SPT with $d L$ =  $d S$ occurs. In this scenario, the laser intensity threshold to induce the SPT is 0.215 Hartree/bohr. Without laser pumping, i.e., $0%$ excited electron, the Born–Oppenheimer molecular dynamic (BOMD) is performed to simulate the structure dynamic at 300 K. When the excited electron is below $0.94%$⁠, the long bond length gradually decreases with the increase of time, while the short bond length steadily increases. At about 130 fs, the long and short bond lengths are closest, thus achieving metastable structures with the strongest metallicity. After that, their difference becomes gradually significant along with the time change. It indicates when the laser intensity is not strong enough to induce the SPT, the crystal transforms from the $M 1$ phase to a metastable structure and then slowly returns to the $M 1$ phase through carrier cooling. These metastable states will be intensively investigated to model THz characteristics of the photoexcited VO $2$ in the following discussions. Similar trends can be also observed from the twisting angle variation [Fig. 2(b)]. Along with the increase of the time, the V–V angles expand to be the largest rapidly and then slowly decrease with the oscillations. Consistently with the changes of bond lengths, the largest angles also increase with the increase of the laser intensity.

When the photoexcited electron exceeds $0.94%$⁠, the evolutions of the bond lengths and the twisting angles with time are illustrated in Figs. 2(c) and 2(d), respectively. The changes of $d L$ and $d S$ above the SPT threshold are different from those below the SPT threshold. Specifically, $d L$ and $d S$ vary to be equal in an ultrafast time scale and then undergo thermal oscillations near $d R$⁠. This physical mechanism is consistent with extensive experiments that the carrier cooling rate after photoexcited SPT is very slow up to a few microseconds.^6,11,40 And the thermal oscillations become increasingly severe as the laser intensity increases. This is mainly because the strong laser injection causes the significant increase of carrier energy. The increase of laser intensity shortens the occurrence moment of SPT from 141 to 56 fs, which are in agreement with the theoretical and experimental time constants, as shown in Table I. The largest twisting angles above the SPT threshold maintain about $176.3 °$⁠, and the subsequent oscillations become more severe as the laser intensity increases, which is consistent with the trends of bond length. The details about the evolution of each V–V bond length and twisting angle with time can be found in Figs. S2 and S3 in the supplementary material.

The designs on dynamically controllable THz devices pay more attention to the metallic state changes of materials. Numerous experiments have verified the separation of IMT and SPT in the photoexcited VO $2$⁠,^6,8,9 with IMT having a faster time scale than SPT, which brings highly attractive ultrafast photoelectric response characteristics. More importantly, the threshold of IMT is lower than that of SPT,^5,11,12 which lays a theoretical foundation for real-time control of the conductivity characteristics of VO $2$⁠.

To investigate these physical mechanisms, we first concentrate on the evolution of the energy level with time to indicate the change of the potential energy surface. Figure 3(a) shows the energy level evolution at the $Γ$ point when $0.73%$ electron is excited. The highest occupied valance band is excited above the Fermi energy at 43 fs, meaning that the crystal transforms to a metal. The photoexcited phonon energy makes the valence energy levels oscillate around the Fermi energy from 43 to 165 fs, maintaining the metallic state. As time increases, the cooling of carrier decreases the photoexcited phonon energy and further results in the gradual decrease of valence energy levels, restoring the insulating state of the crystal. The phenomenon can be demonstrated in terms of density of states (DOS). Figure 3(b) shows the DOS and the projected DOS (PDOS) of the excited carriers at different times when $0.73%$ electron is excited. Obviously, the electrons in the $d ∥$ state are excited into the conduction states by laser pulse. At 40 fs, owing to the highest occupied valance band below the Fermi energy from Fig. 3(a), there is still a bandgap between $d ∥$ states and conduction states, as shown in Fig. 3(b). At 65, 90, and 130 fs within the metallic state interval in Fig. 3(a), $d ∥$ states shift and cross the Fermi energy, combined with the conduction states, and the bandgap disappears accordingly. After 165 fs, the bandgap reopens, as shown in figures at 205 and 209 fs. Similar change trends of the metallic state below the SPT threshold have been reported in experiments.^11,40 It is worth pointing out that the shift of the $d ∥$ state at 130 fs is the largest in Fig. 3(b), which corresponds to the metastable state in Fig. 2(c). We have found that at the moment with the closest $d L$ and $d S$⁠, there are the strongest conductivity characteristic, with the largest shift of the $d ∥$ state. The evolutions of energy level for other excited electron percentages are given in Fig. S4 in the supplementary material, from which two conclusions are drawn, i.e., the IMT moment decreases as the laser intensity increases, and meanwhile, the duration of the metallic state becomes longer. Especially with high intensity photoexcitation (Fig. S4f in the supplementary material with $3.27%$ excited electron), the quasi-instantaneous IMT occurs at 15 fs. The rt-TDDFT + HSE06 hybrid function is further applied to calculate the corresponding time constant, which is also about 15 fs. That further verify the correctness of the rt-TDDFT + $U$ method. The obtained quasi-instantaneous time constants are in agreement with the calculated results¹⁰ and the measured results.⁶ 

Summarizing the time constants for the different excited electrons into Fig. 3(c), the separation with the different thresholds between IMT and SPT can be observed obviously. The IMT emergence is earlier than that of the SPT, with the separation time of 41– 108 fs. The threshold of SPT is $0.94%$ excited electron, while the lower photo excitation (⁠ $0.50%$ excited electron) can induce the IMT. As the laser intensity increases, all time constants including IMT, SPT, and separation time decrease because of the higher photoelectron injection and the higher photoexcited phonon energy, in consistent with the discussions at the beginning of this section.

Figure 3(d) gives the schematic diagram of the $M 1$-to-R phase transition of the photoexcited VO $2$⁠. After the laser pump irradiation, the insulating $M 1$ phase transforms to the metallic state quasi-instantaneously, retaining a weak monoclinic phase. As the laser intensity continuously increases, the metal-like monoclinic phase further transforms to the metallic rutile phase along with time, which signifies the completion of ultrafast SPT.

Unfortunately, the THz spectrum responses of the PIPT system cannot be accurately and efficiently obtained by the traditional DFT. In our previous work,²⁷ a permittivity dispersion modeling scheme of temperature controlled VO $2$ has been proposed, which uses the Drude⁴⁴ model and EMT⁴⁵ to characterize the dispersion properties of the R phase and the intermediate process of IMT. In this work, we extend this method to the photoexcited VO $2$⁠. Since the THz characteristic is sensitive to macroscopic conductivity, the electron transport properties are investigated in this section.

The moments of the strongest metallicity for different laser intensities are listed in Table S2 in the supplementary material. It can be found that below the SPT threshold of $0.94%$ excited electron, the moments of metal-like monoclinic phases are in the range of 129–136 fs. With the increase of the laser intensity above the SPT threshold, the SPT moments decrease from 141 to 56 fs as the previous discussions. The electron transport properties at the moments of metal-like monoclinic phases are solved.

Since the transition of conductivity in the IMT process of VO $2$ is mainly contributed by the change of carrier concentration rather than mobility,^46,47 the carrier relaxation times in the photoexcited phase transition can be approximately regarded as the constant. We use the electron phonon interaction (EPI) method to calculate the relaxation time $τ$ of the R phase at 300 K [Fig. 4(a)], whose calculation details can be found in our previous work.²⁷ The mean $τ=2.68$ fs at the Fermi energy is observed in Fig. 4(a), and this value is applied in the conductivity calculation below.

According to the structure information of the metal-like monoclinic phase including $d L$⁠, $d S$⁠, and $θ$ and SPT moments listed in Table S2 in the supplementary material, we construct the corresponding atomic structures and implement the electron occupations at the corresponding time to solve the Boltzmann transport equations with the calculated relaxation time. The conductivity tensor components for different laser pump intensities are illustrated in Fig. S5 in the supplementary material. The mean conductivities at the Fermi energy are presented in Fig. 4(b) and Table S2 in the supplementary material. The metal-like conductivity variation of photoexcited VO $2$ is approximate linear below the SPT threshold, which is very different from the temperature induced phase transition, where the conductivity undergoes a nonlinear change with an intrinsic hysteresis around $T c$⁠.^2,3,48–50 However, the approximately linear response of conductivity indeed gives rise to the nonlinear variation of the THz characteristic. Hence, the EMT with the nonlinear response is used to characterize the effective dispersion property in the transition process.

The key idea of EMT is using the volume fraction of the metallic phase $f$ to describe the equivalent characteristics in the IMT process. To extend EMT to the PIPT situation, the state of the highest photoexcitation is equivalently considered as the perfect metal state, whose relaxation time $τ$ and conductivity are implemented into the Drude model of the metallic phase. Then, the conductivities of intermediate states in Fig. 4(b) are applied to calculate the corresponding metallic fractions $f$ through the EMT equations. The calculated fractions $f$ for different laser pump intensities below the SPT threshold are illustrated in Fig. 5(a), which vary in a quasilinear way.

To demonstrate the correctness of the proposed ab initio calculations, we further integrate the PIPT EMT into a full-wave electromagnetic solver, i.e., the discontinuous Galerkin time domain (DGTD) method,^51–53 and simulate the THz transmission spectrums of a VO $2$ film sample based on the mica substrate.⁵⁴ The simulation results are compared with the measurements in Fig. 5(b). Great agreement between them indicates the effectiveness of the proposed quantum-electromagnetic computational method. Here, the lower transmission means the stronger metallic state, indicating more parts of the film undergo phase transitions. The measured metallic fractions $f$ are obtained through the inversion of EMT and compared with the theoretical results calculated from ab initio data in Fig. 5(a). The high consistency verifies the correctness of the proposed calculation strategy. Here, the $f$ value closer to 1 indicates there are more parts of the film completing phase transitions. It is worthwhile noting that the used DGTD algorithm is a universal electromagnetic simulation method,^27,51–53 and, thus, the proposed simulation scheme is suitable for arbitrary light-controlled THz device integrated by VO $2$⁠.

In conclusion, this work has explored a quantum-electromagnetic computational method to establish the bridge between micro physical mechanism and macro THz characteristics of photoexcited VO $2$⁠. The dynamic process of quasi-instantaneous IMT and ultrafast SPT is investigated, and the separation with different thresholds of two processes is theoretically discovered, which lays a theoretical foundation for ultrafast photoelectric control. We focus on the electronic transport properties of the metal-like monoclinic phase and use the extended EMT to model the whole PIPT procedure in the THz regime. We believe that the proposed first principles calculation method can not only explain numerous experimental phenomena, but also provide a powerful exploration path and insights for the analysis of the ultrafast mechanism and terahertz properties of other PIPT materials and their optoelectronic THz devices.

See the supplementary material for the bond length evolution, the twisting angle evolution, the evolution of the energy level at Γ point, as well as the volume fractions of the metallic phase for different photoexcited electrons.

We acknowledge that this work was supported in the National Key Research and Development Program of China under Grant No. 2021YFA1401001 and the National Natural Science Foundation of China (NNSFC) (No. 62371355).

The authors have no conflicts to disclose.

Zhen Guo Ban: Data curation (equal); Formal analysis (equal); Investigation (equal); Methodology (equal); Writing – review & editing (lead). Yan Shi: Conceptualization (equal); Data curation (equal); Funding acquisition (lead); Investigation (equal); Methodology (lead); Project administration (lead); Supervision (lead); Writing – original draft (equal); Writing – review & editing (equal). Ning Qian Huang: Formal analysis (equal); Investigation (supporting); Methodology (equal); Writing – review & editing (supporting). Zan Kui Meng: Investigation (supporting); Writing – review & editing (supporting). Shi Chen Zhu: Writing – review & editing (supporting).

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