Over the past three decades it has become clear that the dramatic slowing down dynamics of supercooled liquids is accompanied by a progressive loss of correlation between different kinds of molecular motion. One of the clearest manifestations is the decoupling of translational and rotational diffusion [1, 2]. In high-temperature liquid state, molecular motion is largely homogeneous and follows well-established hydrodynamic laws. The Stokes–Einstein (SE) relation ${D_t} = {k_{\text{B}}}T/C\eta a$ or ${D_{t}}\sim \tau _{{t}}^{ - 1}$ governs the translational diffusion D_t, and the Stokes–Einstein–Debye (SED) relation ${D_r} = {k_{\text{B}}}T/C^{\prime}\eta {a^3}$ or ${D_{{r}}}\sim \tau _{rn}^{ - 1}$ governs the rotational diffusion D_r, where k_B is the Boltzmann constant, T is the temperature, $\eta $ is the bulk viscosity, a is the effective hydrodynamic radius, ${\tau _t}$ is the structural relaxation time, ${\tau _{rn}}$ is the rotational relaxation time, C and C^’ are constants determined by the boundary conditions, where ‘∼’ means ‘proportional’. It implies fundamental coupling between these two diffusion constants. However, the dynamics of liquids become increasingly heterogeneous and non-exponential when undergoing supercooling towards the glass transition point T_g, and the two diffusion constants are decoupled from macroscopic viscosity or relaxation time [3]. The decoupling implies the SE relation or SED relation breakdown, or both fail. By the combination of the two formulas for SE relation and SED relation as well as the assumption of constant a, many different formulas were adopted to verify the decoupling, such as whether ${D_{\text{t}}}/{D_{\text{r}}}$, ${D_{\text{t}}}{\tau _{rn}}$ or ${{T{\tau _{rn}}} \mathord{\left/ {\vphantom {{T{\tau _{rn}}} \eta }} \right. } \eta }$ keep as a constant as the condition changes.

Fujara et al [4] provided the first direct NMR evidence in supercooled ortho-terphenyl (OTP). They observed D_t is proportional to ${\eta ^{ - 1}}$ above 1.2T_g, and otherwise follows a fractional form like ${D_{\text{t}}}\sim {\eta ^{ - 0.75}}$. However, D_r remains proportional to ${\eta ^{ - 1}}$ down to T_g. Mapes et al [5] measured the diffusion constant of OTP near T_g, and found the combination of their data and the data of Fujara can be well fitted by ${D_{\text{t}}}\sim {\eta ^{ - 0.8}}$. Cicerone et al [6, 7] found ${D_{\text{t}}}\sim T/\eta $ and ${D_{\text{t}}}{\tau _r}_n$ are size dependent for the probes in OTP, the SE relation described by ${D_{\text{t}}}\sim T/\eta $ fails with cooling, and ${D_{\text{t}}}{\tau _{rn}}$ changes are relatively small for a large sized probe but increase by almost two orders of magnitude towards ${T_{\text{g}}}$ for the small probe. Eastwood and coworkers [8] reported that decoupling sets in at approximately 1.2T_g based on ${D_{\text{t}}}{\tau _{r2}}$, with quantitative agreement between the simulation and experiment. Stillinger and co-workers [9] observed that the scaled ratio ${D_{\text{t}}}/{D_{\text{r}}}$ is almost equal to 1.0 within 260–346 K, but it started to deviate when T < 260 K. The scaled simulated and experimental ${D_{\text{t}}}/{D_{\text{r}}}$ are almost equal to 1.0 at a high temperature but deviate from 1.0 at a certain temperature upon cooling.

The decoupling phenomena were also intensively investigated in supercooled water. The ${D_{\text{t}}}{\tau _{rn}}$ is observed not to be a constant but is temperature and density dependent in the supercooled SCP/E, TIP5P and ST2 water [10–14]. The decoupling in ${D_{\text{t}}}{\tau _{r2}}$ happens at 2T_g in TIP5P and ST2, but not the usually recognized 1.2T_g, which can be explained by a two-state scenario [15, 16]. The ${D_{\text{r}}}{\tau _t}/T$ is almost a constant in SPC/E water when T > 280 K but starts to increase as the temperature decreases [17]. The ratio ${D_{\text{t}}}/{D_{\text{r}}}$ is not a constant but decreases with cooling in SPC/E [17] and TIP4P/2005 [18]. Both ${D_{\text{r}}}\sim T/{\tau _t}$ and $\tau _{r2}^{ - 1}\sim {T \mathord{\left/ {\vphantom {T {{\tau _t}}}} \right. } {{\tau _t}}}$ get a crossover in ST2 water as the temperature decreases [19]. The $\tau _{r2}^{ - 1}\sim {T \mathord{\left/ {\vphantom {T {{\tau _t}}}} \right. } {{\tau _t}}}$ fails in the whole simulated temperature range; however, ${D_r}\sim {T \mathord{\left/ {\vphantom {T {{\tau _t}}}} \right. } {{\tau _t}}}$ holds at high temperatures and otherwise takes a fractional form.

The similar phenomena are also observed in colloid systems with an increasing volume fraction $\phi $. Weeks and co-workers [20] visualized supercooled colloidal fluids and showed that the D_t began to violate the SE relation prediction at ϕ ≈ 0.52, whereas the D_r did not deviate from the SED relation until ϕ ≈ 0.56–0.57. Subsequent work on quasi-2D colloidal fluids further nuanced this picture, showing that the degree of decoupling is sensitive to local particle geometry, with short dimers exhibiting a different decoupling pattern than long ones, thereby linking the phenomenon to local caging effects [21].

The studies cited above generally report that the translational degree decouples earlier than the rotational one. The prevailing theoretical interpretation attributes this asymmetric decoupling to the spatially heterogeneous and facilitated nature of dynamics in supercooled liquids. In this framework, translation is dominated by rare collective string-like motions that percolate through transient fast regions, allowing it to escape the average viscous drag more readily [22]. In contrast, rotational motion is thought to be more localized and sensitive to the relaxation of the immediate cage formed by neighboring particles, making it more close to the average structural relaxation time and viscosity [23]. This interpretation is strongly supported by molecular dynamics (MD) simulations of simple model liquids, which have visualized these heterogeneous regions and connected the onset of decoupling to distinct stages in the system’s exploration of its potential energy landscape [24].

However, this scenario is observed to be inverted in anisotropic ellipsoids and dumbbells. Quasi-2D colloidal ellipsoids with aspect ratio 6 showed that the ratio D_t/D_r increased by more than one decade across the glass transition, and a lower density of glass transition in the rotational degree is observed than the translational, which implies the rotation decouples earlier than the translation [25, 26]. MD simulations of elongated rods reproduce the same trend when the rotational cooperative length grows faster than the translational one [27]. A dense mixed dumbbells liquid also shows an earlier decoupling of the rotational degree than the translational one at 1.2T_g by comparing ${D_{\text{t}}}\sim \tau _t^{ - 1}$ and ${D_{\text{r}}}\sim \tau _t^{ - 1}$ [22].

These special cases underscore that the order of decoupling is not immutable; rather, it is governed by the relative growth of translational versus rotational dynamic length scales. While the existing experimental or simulated realizations are still sparse, and only exist in anisotropic ellipsoids and dumbbells, a further search for additional model systems is warranted. In this work, we introduce a rigid tetrahedral molecular model as a complementary system to study translation–rotation decoupling. Unlike conventional models such as OTP where the translational diffusion decouples prior to rotational relaxation, our tetrahedral system exhibits an inverse decoupling sequence: rotational dynamics decouple at approximately 2T_g while translational diffusion remains coupled to viscosity and relaxation time to lower temperatures. This inversion originates from the unique isotropic steric hindrance of the tetrahedral geometry, which suppresses rotation more effectively than translation. To our knowledge, this represents the first observation of rotation-preceded decoupling in a tetrahedral molecular glass former.

The configuration for the MD simulation consists of 2048 rigid tetrahedral molecules (TM) in a cubic box with a side length of 9.676 nm. Using the Lennard–Jones potential diameter, and accounting for both the tetrahedral volume and the molecular volume outside the tetrahedron, the volume fraction is approximately $\phi \approx 0.54$. The TM model is based on the coarse-grained OTP model [28, 29]. It consists of four coarse-grained atoms in OTP with the same bond length and interaction parameters as OTP. All our MD simulations were carried out with the GROMACS package [30, 31]. The periodic boundary conditions were applied in all three directions of the Cartesian space. The van der Waals interactions were calculated with a cutoff of 1.4 nm. Thirty four temperatures were simulated and are distributed within 420–1600 K. The temperature was kept at a constant by the Nosé–Hoover thermostat [32, 33].

The D_t is calculated via the mean square displacement as where ${{\boldsymbol{r}}_i}\left( t \right)$ is the center of mass of ith TM at time t, N is the number of TM. The ${\tau _t}$ is determined by the self-intermediate scattering function [34]

where < > denotes a time average, the wavevector $k = 9.0\,{\text{n}}{{\text{m}}^{ - 1}}$ corresponding to the first maximum of the static structure factor. And ${\tau _t}$ is determined by ${F_{\text{s}}}\left( {k,{\tau _t}} \right) = {{\text{e}}^{ - 1}}$. The ${D_{\text{r}}}$ is calculated via its asymptotic relation with the rotational mean square displacement [17, 18] where ${{\boldsymbol{\varphi }}_i}\left( {\Delta t} \right)$ is the angular displacement. The rotational correlation time ${\tau _{rn}}$ is calculated via the rotational correlation function [18, 35] where ${P_n}\left( x \right)$ is the nth order Legendre polynomial, ${{\boldsymbol{e}}_i}\left( t \right)$ is the unit vector of ${{\boldsymbol{\varphi }}_i}\left( t \right)$, and ${\tau _{rn}}$ is determined by ${C_n}\left( {{\tau _{rn}}} \right) = {{\text{e}}^{ - 1}}$.

The $\eta $ is determined by the method proposed by Hess [36] for its reliability and fast convergence in the linear response regime. Since the SE relation is a combination of Einstein relation ${D_{\text{t}}} = {k_{\text{B}}}T/\alpha $ and Stokes’s law $\alpha = C\eta a$. The $\alpha $ is determined by introducing a small force ${f_e}$ on a part of particles in the linear response regime, and 128 TMs are chosen in our simulation. The translational dynamic heterogeneity is usually described by the translational non-Gaussian parameter [37] as

And the translational dynamic heterogeneity is characterized by the rotational non-Gaussian parameter [38] as

The calculated translational diffusion coefficient D_t, rotational diffusion coefficient D_r, viscosity $\eta $, structural relaxation time ${\tau _t}$, rotational relaxation time ${\tau _{rn}}$ for n = 1–6 and frictional coefficient $\alpha $ (scaled by the mass of molecule TM) under different temperature T are plotted in figure 1. The data is usually proposed to follow an Arrhenius law or a VFT law [34]. The fittings show the data can be well described by a VFT law $A = {A_0}{{\text{e}}^{{{ \pm {E_a}} \mathord{\left/ {\vphantom {{ \pm {E_a}} {\left( {T - {T_0}} \right)}}} \right. } {\left( {T - {T_0}} \right)}}}}$ other than the ${D_{\text{r}}}$, which can be well fitted by an Arrhenius law. In practice, we fitted the logarithm of the data to the VFT law in the form $\ln A = \ln {A_0} \pm {{{E_a}} \mathord{\left/ {\vphantom {{{E_a}} {\left( {T - {T_0}} \right)}}} \right. } {\left( {T - {T_0}} \right)}}$. The fitted logarithm of prefactor $\ln {A_0}$, activation energy E_a and the Vogel temperature T₀ are plotted in figure 1 and the fitted data for ${\tau _{rn}}$ are listed in table 1. The fittings show the system behaves more likely as a fragile liquid than a strong one [39]. The trend of D_t, $\eta $ and ${\tau _{r2}}$ versus T are similar to that observed in previous experiments for OTP [4, 5, 40]. McCall et al [40] show that their D_t directly follows the VFT law with E_a = 689 and T₀= 231 K. Their T₀ is close to our value of 220.6 K; however our E_a = 1573.89 is much greater than theirs. The results suggest that the rigid tetrahedral structure either makes the molecular hopping more difficult or requires more free volume for diffusion, consistent with the free volume theory [34]. The T_g is conventionally defined in the literature as the temperature at which the viscosity reaches 10¹² Pa·s [1]. Our VFT fit yields a value of T_g≈ 314.5 K. It is greater than the T_g = 243 K for OTP observed in experiments [5], suggesting that the rigid tetrahedral structure makes the system more prone to glass formation.

To examine whether the translational and the rotational dynamics decouple, we actually need to verify the validity of the SE relation and SED relation. Since the SE relation actually originats from the combination of the Einstein relation ${D_{\text{t}}} = {k_{\text{B}}}T/\alpha $ and Stokes law $\alpha = C\eta a$. If assuming a to be a constant, the SE relation can be described by ${D_{\text{t}}}\sim T/\eta $. And the relation ${D_{\text{t}}}\sim {\tau _t}^{ - 1}$ holds exactly when dynamic heterogeneity is absent. So we adopt the two formulas ${D_{\text{t}}}\sim {\left( {T/\eta } \right)^{ - {\xi _1}}}$ and ${D_{\text{t}}}\sim {\tau _t}^{ - {\xi _2}}$ to test the SE relation, and adopt ${D_{\text{t}}}\sim {\left( {T/\alpha } \right)^{ - {\xi _3}}}$ to test the original Einstein relation. The original SED relation proposed by Debye [35] is ${D_{\text{r}}}\sim \tau _{rn}^{ - 1}$ and ${D_{\text{r}}} = {k_{\text{B}}}T/\varsigma $, the latter can be expressed as ${D_{\text{r}}} = {k_{\text{B}}}T/C^{\prime}\eta {a^3}$ with the Stokes formula $\varsigma = {C^{\prime}}\eta {a^3}$, and is further expressed as ${D_{\text{r}}}\sim T/\eta $, assuming a is a constant. The a is also evaluated by $a\sim {\alpha \mathord{\left/ {\vphantom {\alpha \eta }} \right. } \eta }$ by considering its potential changes. Therefore the SED relation is tested by ${D_{\text{r}}}\sim {\left( {T/\eta } \right)^{{\xi _4}}}$, ${D_{\text{r}}}\sim {\tau _{rn}}^{ - {\xi _{5n}}}$ and ${D_{\text{r}}}{\alpha ^3}\sim {\left( {T{\eta ^2}} \right)^{{\xi _6}}}$. Because the standard errors shown in figure 1 are smaller than 0.1 from 1600 K down to 450 K, and only three data groups within 420–440 K exhibiting slightly larger uncertainties of 0.1–0.2. The fitting is mainly determined by the data in the range 450–1600 K, so if the corresponding exponent ${\xi _i} = 1.0$ (i = 1–6) lies within 0.9–1.1, the SE or SED relation is satisfied within statistical error; otherwise, it is invalid.

As is shown in figure 2, the SE relation ${D_{\text{t}}}\sim T/\eta $, ${D_{\text{t}}}\sim \tau _t^{ - 1}$ and ${D_{\text{t}}}\sim T/\alpha $ follow fractional forms as ${D_{\text{t}}}\sim {\left( {T/\eta } \right)^{{\xi _1}}}$, ${D_{\text{t}}}\sim {\tau _t}^{ - {\xi _2}}$ and ${D_{\text{t}}}\sim {\left( {T/\alpha } \right)^{{\xi _3}}}$; and the exponents are ${\xi _1} = 0.956$, ${\xi _2} = 1.029$ and ${\xi _3} = 1.023$, respectively. All exponents are so close to the exact result ${\xi _i} = 1.0$, which indicates the SE relation is valid in the simulated temperature range. The ${D_t}\sim \tau _t^{ - 1}$ is an exact result when the translation displacement $\Delta {{\boldsymbol{r}}_i}\left( t \right)$ follows Gaussian distribution. However, Kawasaki and Kim have shown that the validity of ${D_{\text{t}}}\sim \tau _t^{ - 1}$ depends on the ratio of the cage-breaking time scale to structural relaxation [41]. The cage-breaking time is characterized by t_max, defined as the time corresponding to the maximum of the ${\alpha _{2t}}\left( t \right)$. When the ratio ${{{t_{\max }}} \mathord{\left/ {\vphantom {{{t_{\max }}} {{\tau _t}}}} \right. } {{\tau _t}}}$ is not much smaller than unity (e.g. ∼0.1), diffusion and relaxation remain coupled and the ${D_{\text{t}}}\sim \tau _t^{ - 1}$ is preserved. As shown in figure 3(c), the smallest ratio at T = 420 K is so close to 0.1. This explains why, despite the maximum ${\alpha _{2t}}\left( t \right)$ being around 1.5 at T = 420 K as shown in figure 3(a), the ${D_{\text{t}}}\sim \tau _t^{ - 1}$ remains valid. The discrepancies in exponents between ${D_{\text{t}}}\sim {\left( {T/\eta } \right)^{{\xi _1}}}$ and ${D_{\text{t}}}\sim {\left( {T/\alpha } \right)^{{\xi _3}}}$ are due to the assumption of constant a adopted in ${D_{\text{t}}}\sim T/\eta $. By comparing ${D_{\text{t}}}\sim {\left( {T/\eta } \right)^{{\xi _1}}}$ and ${D_{\text{t}}}\sim {\left( {T/\alpha } \right)^{{\xi _3}}}$, the a should be varied with D_t like $a\sim D_{\text{t}}^{0.0685}$ and decreases with decreasing temperature. In case of the small differences between ${\xi _1}$ and ${\xi _3}$, the temperature dependence is very small. Nevertheless, these modest variations still induce discrepancies of approximately 0.08 between the exponent ${\xi _1}$ and ${\xi _3}$, and which propagate into much larger differences in the subsequent test of the SED relation. This suggests that exponent a serves as a sensitive diagnostic for the validity of the SE relation. Such behavior aligns with our earlier observation that a exhibits a downward trend upon cooling in systems dominated by attractive interactions [42].

The simulated results of SED relation ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ plotted in figure 4(a) do not fall on a line but with a crossover at temperature ${T_x} \approx 600\,{\text{K}}$. The exponents for the two parts are ${\xi _4} = 0.818$ and 0.289, respectively. Both fitted exponents depart from the theoretical value of ${\xi _4} = 1.0$, lying beyond the validity interval 0.9–1.1. Notably, the deviation is significantly attenuated in the high-temperature regime T> T_x compared with the low-temperature regime T < T_x. The result indicates the breakdown of the SED relation described by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ in the whole simulated temperature range.

A similar crossover is also observed around 600 K in the SED relation described by ${D_r}\sim \tau _{rn}^{ - 1}$ for n = 1–6, but with different exponents compared to ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$. As the fitted ${\xi _{5n}}$ listed in table 2 is shown, the exponents in ${D_r}\sim {\tau _{rn}}^{ - {\xi _{5n}}}$ for n = 1–5 is so close to the exact result ${\xi _{5n}} = 1.0$ above 600 K, which indicate the validity of ${D_r}\sim \tau _{rn}^{ - 1}$. However, ${\xi _{5n}} \approx 0.6$ implies the breakdown of ${D_r}\sim \tau _{rn}^{ - 1}$ for n = 1–5 below 600 K. Moreover, the ${D_r}\sim \tau _{r6}^{ - 1}$ is breakdown in the whole temperature range and with ${\xi _{56}} = 1.227$ for T > 600 K and otherwise 0.585. The result is consistent with ${\alpha _{2r}}\left( t \right)$ plotted in figure 3(b), which is close to Gaussian distribution at higher temperatures but deviates as temperature decreases. The results are similar with observations in TIP4P/2005, where the ${D_r}\sim \tau _{rn}^{ - 1}$ is more strongly violated at larger n [18]. As pointed out in [18], lower order n probes large angle molecular reorientation, while high-order probes small-angle motion within the cage, the two are fundamentally different physical processes with distinct time scales. At high temperatures, molecules rotate freely, cage effects are minimal, and both processes remain coupled to thermal motion. Lower order ${\tau _{r1 - 5}}$ follows Debye relaxation, matching the ${D_r}$, so ${D_r}\sim \tau _{r1 - 5}^{ - 1}$ holds. At low temperatures, dynamic heterogeneity emerges. Cage escape becomes rare and heterogeneous, slowing ${\tau _{r1 - 5}}$ dramatically. However, n = 6 still probes fast in-cage wobbling. Thus ${\tau _{r1 - 5}}$ and ${D_r}$ decouple, ${D_r}\sim \tau _{r1 - 5}^{ - 1}$ breaks down, while ${\tau _{r6}}$ was never synchronized with ${D_r}$ due to this inherent time-scale separation.

The ${\alpha _{2r}}\left( t \right)$ is larger than ${\alpha _{2t}}\left( t \right)$ at any temperature, which indicates the system has more heterogeneous rotational dynamics than the translational dynamics at any temperature. This arises from steric hindrance, as observed in colloidal ellipsoids [25, 26] and elongated rods [27]. The rigid tetrahedral structure and high volume fraction hinder molecular rotation more severely than translation, since steric effects require neighboring molecules to rotate collectively to accommodate the central molecule’s reorientation.

The data for the SED relation ${D_r} = {{{k_{\text{B}}}T} \mathord{\left/ {\vphantom {{{k_{\text{B}}}T} {{C^{\prime}}\eta {a^3}}}} \right. } {{C^{\prime}}\eta {a^3}}}$ tested by ${D_r}{\alpha ^3}\sim {\left( {T{\eta ^2}} \right)^{{\xi _6}}}$ plotted in figure 4(c) is also shown as a crossover around ${T_x} \approx 600{\text{K}}$. The fitted exponent is ${\xi _6} = 1.051$ above T_x and otherwise ${\xi _6} = 1.262$. The former is close to the exact ${\xi _6} = 1.0$ and indicates the ${D_r} = {{{k_{\text{B}}}T} \mathord{\left/ {\vphantom {{{k_{\text{B}}}T} {{C^{\prime}}\eta {a^3}}}} \right. } {{C^{\prime}}\eta {a^3}}}$ is established when T > T_x. The latter lies beyond the validity interval and implies the breakdown for T < T_x. Comparing the results given by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ and ${D_r} = {{{k_{\text{B}}}T} \mathord{\left/ {\vphantom {{{k_{\text{B}}}T} {{C^{\prime}}\eta {a^3}}}} \right. } {{C^{\prime}}\eta {a^3}}}$, the results given by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ take some corrections after considering the variation of a into account. The ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ is almost corrected to the exact result ${\xi _6} = 1.0$ for T > T_x. However, the ${\xi _4} = 0.289$ for T< T_x is corrected to ${\xi _6} = 1.262$. Although the correction makes the breakdown smaller, yet the ${D_r} = {{{k_{\text{B}}}T} \mathord{\left/ {\vphantom {{{k_{\text{B}}}T} {{C^{\prime}}\eta {a^3}}}} \right. } {{C^{\prime}}\eta {a^3}}}$ is still invalid for T< T_x. The result also signifies the temperature dependent of a and its significance in testing of the validity SED relation.

Combining the results given by figures 2 and 4, the SE relations are all valid for the three forms but the three formulas of the SED relation are all breakdown. The decoupling is existed in the whole simulated temperature range by comparing the three SE relations with the SED relation described by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ and ${D_r}\sim \tau _{r6}^{ - 1}$. Comparing the results given by the three SE relations with the SED relation described by ${D_r}\sim \tau _{rn}^{ - 1}$ for n = 1–5 and ${D_r} = {{{k_{\text{B}}}T} \mathord{\left/ {\vphantom {{{k_{\text{B}}}T} {{C^{\prime}}\eta {a^3}}}} \right. } {{C^{\prime}}\eta {a^3}}}$, the translational and rotational motion are still coupling for T> T_x but the two decouple when T< T_x. The T_x and T_g nearly satisfy the relation ${T_x} \approx 2{T_{\text{g}}}$, which is similar to that observed in TIP5P and ST2 water [13]. However, it is different from that observed in OTP [4, 9] and dumbbell molecular systems [22], which shows a decoupling at 1.2T_g. Our systems show a similar earlier rotational decoupling than the translational motion as that observed in dumbbell molecular systems, but inverted in OTP [4, 9], TIP5P and ST2 water [13].

The phenomena observed can be explained by the differences in the mechanism of diffusion between the translational and the rotational [22, 23]. Translational diffusion is nonlocal, requiring escape from the cage trap. In contrast, rotational motion is more localized and sensitive to the relaxation of the immediate cage formed by neighboring particles. As observed by Kob and coworkers [22, 43] for a two-component dumbbell mixture at $\phi = 0.708$, rotational decoupling occurs earlier in this system despite its weak steric hindrance arising from the dumbbell geometry. This system resembles commonly simulated water or OTP, differing primarily in its significantly higher volume fraction. When steric hindrance and dense packing are present, rotation is typically more constrained than translation, leading to earlier decoupling of the rotational motion. Our system has a large volume faction $\phi = 0.54$, and the tetrahedral structure introduces strong steric hindrance to rotation. We can make a simple estimate of the steric hindrance for rotation. Assuming rotation about the center of mass, both the tetrahedral side length and the Lennard–Jones parameter σ are equal to 0.483 nm, each molecule requires a circumscribed sphere of a tetrahedron of 0.65nm³ to rotate freely, and the total volume need is approximately 1.47 times the simulation box volume. In comparison, for translation, the molecule has an additional 0.46 fraction of free volume available for motion. Our rough estimate shows that steric hindrance has a stronger effect on rotation than on translation, leading to earlier rotational decoupling.

To further explore the possible influence of steric hindrance introduced dynamic heterogeneity on the translational and rotational dynamics, we calculated the D_t, D_r, ${\tau _t}$ and ${\tau _{r3}}$ for the fastest 7% and slowest 7% of molecules, following the approach in [17], and then examined the corresponding SE relation ${D_t}\sim \tau _t^{ - 1}$ and SED relation ${D_r}\sim \tau _{r3}^{ - 1}$, the correlated data are plotted in Fig. 5. The D_t and ${\tau _t}$ for the fastest and slowest subsets can differ by up to 50 fold, indicating strongly heterogeneous translational dynamics. The ${D_t}\sim \tau _t^{ - 1}$ holds for the slowest 7% of molecules as shown by ${\xi _s} = 1.004$ but breaks down for the fastest with ${\xi _f} = 1.93$. While it has been argued that D_t is determined by mobile molecules whereas ${\tau _t}$ is governed by immobile ones [34], we examined a mixed SE relation ${D_t}\sim \tau _t^{ - 1}$ using the D_t of the fastest and the ${\tau _t}$ of the slowest. Despite the large discrepancy between these two subsets, the SE relation remains valid in this mixed case as shown by the green line in figure 5(a), consistent with the results shown in figure 2.

Compared with the D_t and ${\tau _t}$ for the fastest and slowest subsets, no such large difference is observed in D_r or ${\tau _{r3}}$, and these differ by only up to 1.8 fold, suggesting that the rotational dynamics are more homogeneous than the translational ones. However, the SED relation ${D_r}\sim \tau _{r3}^{ - 1}$ for both the fastest and slowest subsets is only valid at high temperatures with ${\xi _f} = 0.925$ and ${\xi _s} = 1.02$, and breaks down at low temperatures as shown by ${\xi _f} = 0.853$ and ${\xi _s} = 0.854$, consistent qualitatively with the results shown in figure 4. The mixed SED relation ${D_r}\sim \tau _{r3}^{ - 1}$ using the D_r of the fastest and the ${\tau _{r3}}$ of the slowest fails in the whole simulated temperature range with ${\xi _{fs}} = 0.871$. Given that the validity of the ${D_t}\sim \tau _t^{ - 1}$ and ${D_r}\sim \tau _{r3}^{ - 1}$ relies on Gaussian behavior, steric hindrance causes greater deviation from Gaussian dynamics in rotation than in translation, consistent with the observed non-Gaussian parameter.

The $\phi = 0.54$ is close to the decoupling displayed in the spherical colloid system [20]. However, due to the rigid tetrahedral structure discussed above, each molecule requires a large spherical volume to rotate freely. In contrast, the decoupling in anisotropic particles such as dumbbells, ellipsoids, and elongated rods is due to the intrinsic translation–rotation coupling, prolate molecules can translate along their long axis via rotation, leading to a shape-dependent decoupling that varies with aspect ratio and packing fraction. Thereby it requires less free volume for rotational motion than our system. Consequently, these systems need a higher volume fraction to achieve earlier rotational decoupling relative to translation. By comparison, the effective volume fraction in our system is significantly lower than that in dumbbell systems [22, 43], colloidal ellipsoids [25, 26] and elongated rods [27], where the volume fraction $\phi > 0.7$. The steric hindrance with rigid tetrahedral structure is also the reason for the decoupling exhibiting temperature differences between our systems and the dumbbell system [22, 43]. Meanwhile, the decoupling temperature of approximately 2T_g is close to that observed in TIP5P and ST2 water models [13, 15, 16], which might suggest that the decoupling in TM is similarly related to its tetrahedral molecular structure. However, the underlying mechanisms are fundamentally different. In supercooled water, hydrogen bonding drives molecules into tetrahedral configurations upon cooling, triggering a disorder-to-order transition between two distinct liquid states. The emergence of a system–spanning hydrogen-bonded network is captured by the tetrahedral order parameter Q₄, which exhibits a pronounced increase as the transition proceeds. By contrast, our system features only isolated TM.

To uncover the hidden resemblance, we calculated the local Q₄ based on the molecular center of mass [44]. The corresponding probability distribution and ensemble-average mean are present in figure 6. The probability distributions differ only slightly across the entire temperature range, from 1600 K down to 420 K. The mean Q₄ exhibits only a weak temperature dependence, increasing by merely 0.03 upon cooling from 1600 K to 420 K. Notably, no liquid–liquid phase transition is detected, in contrast to the behavior reported for supercooled water [44]. Moreover, previous studies have shown that the decoupling temperature is not universal but depends on the aspect ratio of rods [45], as well as the bond length of dimers and the volume fraction [21]. Therefore, we propose the similar decoupling temperatures observed in our system and supercooled water may be coincidental.

In summary, we have examined the coupling and decoupling between the translational and rotational dynamics in a tetrahedral molecular liquid by testing the validity of the SE relation and SED relation. Our results indicate the SE relation described by ${D_t}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$, ${D_t}\sim \tau _t^{ - 1}$ and ${D_t}\sim {T \mathord{\left/ {\vphantom {T \alpha }} \right. } \alpha }$ are all established. The SED relation described by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$, ${D_r}\sim \tau _{rn}^{ - 1}$ and ${D_r}\sim {T \mathord{\left/ {\vphantom {T {\eta {a^3}}}} \right. } {\eta {a^3}}}$ all show a crossover around T_x = 600 K. The ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ is breakdown in the whole simulated temperature range and follows a fractional form on both sides of T_x. The ${D_r}\sim \tau _{rn}^{ - 1}$ for n = 1–5 are all valid for T> T_x but is invalid for ${D_r}\sim \tau _{r6}^{ - 1}$. The ${D_r}\sim \tau _{rn}^{ - 1}$ for n= 1–6 are all breakdown when T < T_x and in a factional form with ${\xi _{5n}} \simeq 0.6$. The ${D_r}\sim {T \mathord{\left/ {\vphantom {T {\eta {a^3}}}} \right. } {\eta {a^3}}}$ tested by ${D_r}{\alpha ^3}\sim {\left( {T{\eta ^2}} \right)^{{\xi _6}}}$ is also established when T> T_x and otherwise breakdown. The differences between ${\xi _1}$ and ${\xi _3}$ imply the a is temperature dependent, although the dependence is small. However, the result given by ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ is almost corrected to the exact result ${\xi _6} = 1.0$ when T> T_x and the breakdown of ${D_r}\sim {T \mathord{\left/ {\vphantom {T \eta }} \right. } \eta }$ is also reduced after taking the variation of a into account. It implies the importance of considering the variation of a in testing the SE relation and SED relation. Combining the results given by the SE relation and SED relation, our simulations indicate the translational motion is still coupling with the rotational motion for T> T_x and otherwise decouple. The results are due to the steric hindrance introduced by the rigid tetrahedral structure and dense packing, which makes the rotation more strongly hindered than the translation. Overall, the central finding of this study is the identification of an inverse hierarchy of dynamic arrest in a rigid tetrahedral model: contrary to OTP and most existing glass formers where translation decouples first, rotation decouples earlier than translation in our system. This demonstrates that molecular geometry can fundamentally alter the competition between translational and rotational degrees of freedom in supercooled liquids. However, there are still some questions that need to be answered in the future; as we have discussed above, the volume fraction and steric hindrance have a strong influence on the decoupling order, whether there exists a critical volume fraction at which the earlier decoupling of translational dynamics relative to rotational dynamics will be reversed, like in supercooled water or OTP, especially for spherical particle system.