Quantum structural fluxion in superconducting lanthanum polyhydride

Article https://doi.org/10.1038/s41467-023-37295-1 Quantum structural fluxion in superconducting lanthanum polyhydride Hui Wang 1,2 ,PascalT.Salzbrenner 3 , Ion Errea 4,5,6 , Feng Peng 7 , Ziheng Lu 3 ,HanyuLiu 2,8 , Li Zhu 9 ,ChrisJ.Pickard 3,10 &YansunYao 11 The discovery of 250-kelvin superconducting lanthanum polyhydride under high pressure marked a significant advance toward the realization of a room‐ temperature superconductor. X-ray diffraction (XRD) studies reveal a nonstoichiometric LaH 9.6 or LaH 10±δ polyhydride responsible for the superconductivity, which in the literature is commonly treated as LaH 10 without accounting for stoichiometric defects. Here, we discover significant nuclear quantum effects (NQE) in this polyhydride, and demonstrate that a minor amount of stoichiometric defects will cause quantum proton diffusion in the otherwise rigid lanthanum lattice in the ground state. The diffusion coefficient reaches ~10−7cm2/s in LaH 9.63 at 150 gigapascals and 240 kelvin, approaching the upper bound value of interstitial hydrides at comparable temperatures. A puzzling phenomenon observed in previous experiments, the positive pressure dependence of the superconducting critical temperature T c below 150 gigapascals, is explained by a modulation of the electronic structure due to a prematuredistortionofthehydrogen lattice in this quantum fluxional structure upon decompression, and resulting changes of the electron-phonon coupling. This finding suggests the coexistence of the quantum proton fluxion and hydrogen-induced superconductivity in this lanthanum polyhydride, and leads to an understanding of the structural nature and superconductivity of nonstoichiomectric hydrogen-rich materials. Superconductivity at near room temperature has been discovered in clathrate polyhydrides at megabar pressures1–7. Determining the crystal structure responsible for the superconductivity is of critical importance, yet a great challenge8,9. Experimentally, difficulties in probing lighter elements arise in XRD, impeding a direct determination of the complete lattice symmetry. The measured T c and its pressure dependence (dT c /dp) tighten constraints on structural models; however, calculations on the candidate structures based on Bardeen–Cooper–Schrieffer (BCS) theory10 unequivocally suggest a negative dT c /dp11–20, which contradicts the experimentally observed Received: 29 November 2022 Accepted: 9 March 2023 Check for updates 1 Key Laboratory for Photonic and Electronic Bandgap Materials (Ministry of Education), School of Physics and Electronic Engineering, Harbin Normal University, 150025 Harbin, China. 2 International Center for Computational Method & Software, College of Physics, Jilin University, 130012 Changchun, China. 3 Department of Materials Science & Metallurgy, University of Cambridge, 27 Charles Babbage Road, Cambridge CB3 0FS, UK. 4 Fisika Aplikatua Saila, Gipuzkoako Ingeniaritza Eskola, University of the Basque Country (UPV/EHU), Europa Plaza 1, 20018 Donostia/San Sebastián, Spain. 5 Centro de Física de Materiales (CSIC-UPV/EHU), Manuel de Lardizabal Pasealekua 5, 20018 Donostia/San Sebastián, Spain. 6 Donostia International Physics Center (DIPC), Manuel de Lardizabal Pasealekua 4, 20018 Donostia/San Sebastián, Spain. 7 College of Physics and Electronic Information, Luoyang Normal University, 471022 Luoyang, P. R. China. 8 State Key Laboratory of Superhard Materials and International Center of Future Science, Jilin University, 130012 Changchun, China. 9 Department of Physics, Rutgers University, Newark, NJ 07102, USA. 10 Advanced Institute for Materials Research, Tohoku University 2-1-1 Katahira, Aoba, Sendai 980-8577, Japan. 11 Department of Physics and Engineering Physics, University of Saskatchewan, Saskatoon, Saskatchewan S7N 5E2, Canada. e-mail: [email protected] Nature Communications | (2023) 14:1674 1 1234567890():,; 1234567890():,; positive dT c /dpin some pressure regimes (Supplementary Fig. 1). This phenomenon was initially observed in the 250-kelvin superconducting lanthanum polyhydride1, and later in polyhydrides of yttrium6,7and calcium21 with T c reaching 257 K and 212 K, respectively. This ‘positive dT c /dpcontradiction’presents a major obstacle to our full understanding of the crystal structures of these superconducting polyhydrides. As the first superconductor with a T c above 250K, lanthanum polyhydride has been studied by several groups independently1–5.In addition to the high T c , a positive dT c /dpwas observed between 137 and 150GPa by measurements on four samples (samples 1–4inRef.1) synthesized at different pressures with laser heating. XRD measurementsdeterminethatthe lanthanumatomsforma face-centered cubic (fcc) lattice at pressures of 137–218 GPa, while the hydrogen atoms haveundeterminedlocationswithinthelattice.Basedonthemeasured crystal volume, the hydrogen-to-lanthanum (H/La) ratio around the maximum T c was estimated to be 9.6 (150GPa)1or 9–11 (180–200 GPa)2,3in the two studies, respectively. This trend was recently confirmed by new experiments on the fcc →C2/mphase transformation at p c = 135 GPa, which shows an even steeper decrease of T c below p c 5. With a moderate synthetic pressure and high symmetry, ‘fcc’ lanthanum polyhydride provides a model system for exploring the structural nature of high-T c hydrides in both experiment and theory. In fact, ab initio calculations have guided the experimental discovery of lanthanum polyhydrides, and predicted the appearance of high T c superconductivity, i.e., an estimated T c of 280-kelvin for fcc-LaH 10 at 210 GPa11,12. This prediction turned out to be very close to the fcc lanthanum polyhydride later synthesized, with some differences in T c , synthetic pressure, and hydrogen content. Recently, it has been theorized that the inclusion of quantum atomic fluctuation is essential for a correct calculation of T c and the pressure boundary of fcc-LaH 10 , which is dictated by the quantum nature of lanthanum polyhydride structures16.However,anegativedT c /dppredicted at 137–150 GPa16 remains in contradiction to the experimental observations (Supplementary Fig. 1). This disagreement suggests that critical factors have been overlooked in previous calculations of superconductivity. In particular, the stoichiometric defect observed by experimental studies may play an important role in determining the structure and underlying superconductivity. Results and discussion Although minor in amount, defects in solid materials can strongly affect their properties. The thermodynamic stability of a defect in high-pressure solids as well as the relative stabilities of different defect structures can be evaluated by the defect formation enthalpy (Hf). As shown in Fig. 1a, the crystal structure of fcc-LaH 10 is represented as the insertion of an ‘Hcube’into octahedral interstices of the fluorite-LaH 2 structure. In this model, we calculated the Hffor a vacancy defect either at a corner of the H cube (V C ) or at a tetrahedral interstice of the La lattice (V T ), and found that the configurationally averaged Hfbecomes negative below 158 GPa. This suggests that fccLaH 10 is prone to vacancy defects, which agrees well with the experimental finding of hypostoichiometric LaH 10-δ at 150 GPa specifically LaH 9.6 1. Notably, the calculated ‘vacancy occurring region’ coincides the region where the T c of the fcc lanthanum polyhydride depends positively on the pressure (which previous calculations have failed to predict), indicating that the vacancy structure has a key role to play in the superconductivity. A low concentration of vacancies in hydrogen sites does not change the XRD pattern of the fcc lanthanum polyhydride, which is determined primarily by the La sublattice. However, these vacancies significantly affect the dynamics of hydrogen in the crystal (Supplementary Fig. 2). The distortion in the hydrogen sublattice spreads out away from the vacancy sites in both the V C and V T models of LaH 10 with an H/La ratio of 9.97, and even results in a ‘liquid-like’Hframeworkin 135 140 145 150 155 160 165 -60 -50 -40 -30 -20 -10 0 10 20 135 150 165 180 31 32 33 34 35 b VT VC Average Hf (meV) Pressure (GPa) 158 GPa a Expt. Drozdov LaH9.6 LaD10 Cal. This work LaH10 LaH9 LaD10 VCMD Volume ( /f.u.) Pressure (GPa) VT VC H La Å3 Fig. 1 | Vacancy formation enthalpy and pressure-volume relation. a The formation enthalpy (H f ) of a single vacancy in fcc-LaH 10 at two inequivalent lattice sites in the clathrate hydrogen framework, V T and V C ,andtheirconfigurational average calculated at different pressures using 2 × 2 × 2 extension of conventional unit cell of fcc-LaH 10 . The locations of V T and V C are illustrated in a conventional unit cell of the fcc-LaH 10 structure by red and black balls, respectively. bThe experimental pressure-volume relations of fcc-LaH 9.6 and fcc-LaD 10 measured by Drozdov et al.1 compared to the theoretical values of fcc-LaH 10 ,fcc-LaH 9 and fcc-LaD 10 derived from quantum simulations at constant pressure and temperature of 300K. The volumes selected for subsequent quantum simulation (V CMD )forfcc-LaH 9.63 are thereby linked to the ‘quantum’pressure. Formula unit is abbreviated as f.u. Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 2 the former at 150GPa. For LaH 10-δ , the presence of vacancies changes the potential energy surface, while quantum nuclear fluctuation is expected to affect the transport behavior of vacancies. The calculated zero-point energy (ZPE) of ~150 meV/LaH 10 is comparable to the vacancy migration energy of 160/200 meV for the migration path of V C →V C’ /V C →V T at 150GPa. This suggests that nuclear zero-point motion can promote proton hopping to neighboring vacancies and cause the LaH 10-δ ground state to not be a single structure around which the atoms vibrate, but instead one where the protons dynamically explore different vacancies in a fixed fcc-La framework. We refer to this behavior as being ‘fluxional’in the sense used by Goncharov et al. for phase IV of solid hydrogen22. It is well-established by experiment and theory that thermal effects dominate the fluxionality in this structure: A classically ordered phase III replaces the fluxional phase IV at low temperature23–25. Here, we demonstrate that, in the presence of vacancies, hydrogen NQEs are strong enough—and indeed essential—to stabilize a fluxional framework in the fcc lanthanum polyhydride. Therefore, we are dealing with a quantum fluxional structure (QFS). Taking LaH 9.63 as an example, we investigate the NQE and thermal effects on the structuralfluxion at 150 GPa and low temperatures up to 240 K. The stoichiometry is selected according to the experimental estimation (LaH 9.6 ). Our theoretical simulations suggest a H/La ratio of 9.54 (or 9.71) for the experimental samples synthesized at 150 GPa (or 137 GPa), in good agreement with experiments1(Fig. 1b). We adopt ab initio centroid molecular dynamics (CMD)26 to treat the nuclei quantum mechanically. As shown by the mean square displacement (MSD) curves in Fig. 2a, the simulations reveal appreciable proton migration, with equivalent mean migration distances of vacancies reaching 1.0 to 2.0 Å between 60 and 240 K (Supplementary Fig. 3). These distances arealmost aslargeasorlargerthanthe shortestH-Hseparation around 1.2Å.Thisisconsistentwiththe‘network’-like density distribution patterns obtained, for instance, at 60 K (Fig. 2b). In contrast, the diffusivity is much smaller in ab initio molecular dynamics (MD), where the nuclei are treated classically (Fig. 2a, c). This confirms that the structural fluxion in LaH 9.63 is dominated by the NQE, rather than thermal effects. The crucial role of vacancies in this intriguing quantum phenomenon is demonstrated by a comparison to the dynamics of stoichiometric fcc-LaH 10 at the same pressure, which exhibits no diffusion below 700 K27. Experimentally, the occurrence of ‘quenched-in’vacancies in polyhydrides is inevitable, as a result of the thermal treatment of the samples. This finding therefore suggests that vacancy effects interact with the strong NQE to result in strong quantum structural fluxion in polyhydrides. Weapproximatethe proton diffusion coefficient DinLaH 9.63 from the slope of the MSD in simulations at 240K and between 137 to 176 GPa. The temperature is chosen to be close to the measured T c of the fcc lanthanum polyhydride. As shown in Fig. 3a, the Dvalue falls in the range of 10−6to 10−7cm2/s in our quantum simulation, whereas it is significantlylower in the classical simulation. The Dvalue is 2–3orders of magnitude below the threshold criterion for classical superionicity (~10−4cm2/s) for freely diffusing protons28. This notwithstanding, it reaches or indeed exceeds the upper bound diffusivities observed in interstitial hydrides at room temperature, such as 3.8 × 10−7cm2/s in fcc-Cu 2 H29, and approaches that of phase IV hydrogen30.With increasing pressure, Ddecreases in interstitial hydrides (e.g. Cu 2 H29 and FeH31), due to the contraction of the metallattice and the resulting increase of the activation energy for proton hopping. In phase IV hydrogen, on the other hand, Dincreases with pressure30, probably owing to the drastic increase of the proton hopping rate with the shortened of the distance. Our preliminary results suggest that the unique ‘host-guest’structure of LaH 9.63 likely induces a significant competition between the two mechanisms mentioned, resulting in a nonmonotonic pressure trend of the coefficient D. In LaH 9.63 , proton diffusion breaks the balance of internal Coulomb repulsions in the ‘Hcube’located in octahedral interstices of La sublattice, which therefore distorts the H sublattice. Lattice expansion promotes the distortion, as indicated by the heightening of the second coordination shell (peaked at ~1.85 Å) in the average radial distribution function [RDF, g(r), Supplementary Fig. 4], which is correlated to a smearing of the interand inner-cube H-H separations. Based on a configurational distance (denoted as ξ)ofcrystalfingerprint matrices32, which is sensitive to local structural changes, we parametrized the structural difference of the QFS relative to the crystal lattices of static fcc and C2/mphases, and a triclinic P1 structure mimicking the QFS, as showninFig.3bbyξ qf ,ξ qc and ξ qp respectively. For the H substructure, we find ξ qf >ξ qc >ξ qp ,withtheξ qf much beyond the other two at larger volumes, revealing an increasing distortion of the H sublattice in QFS from fcc to lower symmetry upon decompression. For the La substructure, ξ qc >ξ qf ≈ξ qp ,withξ qc much larger than ξ qf and ξ qp .This agrees with the XRD measurements that determine an fcc sublattice for La atoms above 137 GPa1. The tiny distortions related to ξ qf and ξ qp are both within the uncertainty of refinements for the fcc phase in XRD studies (Supplementary Fig. 4). The results suggest a premature lowsymmetry distortion of the H substructure upon volume expansion relative to the La substructure in the quantum fluxional LaH 9.63 . Using the centroid configurations of CMD simulations at 240 K, i.e. the QFS, we calculated the electronic density of states at the Fermi level NϵF  in LaH 9.63 , which exhibits opposite pressure trends at pressure below and above 150 GPa (Fig. 4a). The fcc-LaH 10 (quantum) crystal has a monotonously negative pressure dependence of NϵF  above 100 GPa12,16, and this prediction can be extended to LaH 9.63 using the rigidband model of the electronic structure by virtue of an artificial shift of ϵF33. The pressure dependences of NϵF  in the C2/mphase of LaH 10 and LaH 9.63 (within the rigid-band approximation) are both non-monotonic, similar to quantum fluxional LaH 9.63 , suggesting that pressure effects on dNϵF  /dparemuchmoresignificant in distorted structures compared to the high-symmetry fcc phase. The non-monotonic pressure trend of NϵF  in LaH 9.63 is a statistical consequence of the QFS having an average fcc-La sublattice with diversely distorted H substructures. This trend cannot be attributed to a single classical configuration, but can roughly be reproduced by the P1 structure mimicking the QFS 01234 0.00 0.04 0.08 0.12 0.16 0.20 0.24 Classical Quantum Nuclear density distributions Quantum: 240 K 120 K 60 K Classical: 240 K 120 K 60 K MSD ( ) Time (ps) Å2 ab c Fig. 2 | Structural fluxion in LaH 9.63 at 150 GPa. a The proton MSD derived from centroid trajectories of the CMD simulations and those of MD simulations. bThe [100] view of the quantum nuclear density distribution at 60 K extracted from a CMD simulation, with full consideration of the 16 beads. Neighboring protons are illustrated in different colors. cthe [100] view of the classical nuclear density distribution in 16 MD simulation runs of 4 ps distinguished by various proton colors. The 16 runs were initialized from different centroid configurations of the 4-picosecond CMD trajectory with a sampling interval of 0.25 ps. Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 3 (Supplementary Fig. 5), which implies that the premature distortion of the H substructure (e.g. fcc →P1 symmetry) accounts for the sign change of dNϵF  /dpin the quantum fluxional LaH 9.63 . The analytic McMillan34 and Allen-Dynes (AD)35 formulas (see details in Supplemental Information) have played an important role in analyzing the mechanism of pressure-dependent superconductivity in high-T c hydrides. In these formulas, the electron-phonon coupling (EPC) constant λdepends explicitly on characteristic parameters of both electrons and phonons, in addition to the average of the electron-phonon matrix elements, hI2i.Specifically, λcan be expressed as the product of hI2i=M (denoted as β,withMbeingatomicmass)andNϵF  = ω2 2  (denoted as ζ,with ω2being the second frequency moment of the Eliashberg function). ζdirectly describes the competition between electron and phonon contributions. Analysis of literature data suggest that ζdecreases more than two times faster than βincreases upon compression (Supplementary Table 1). Furthermore, it was found that dλ/dpplays a dominating role in determining the dT c /dp,despitethefactthattheT c also depends explicitly on ωlog and  ω2/ωlog (with ωlog being the logarithmic frequency moment of the Eliashberg function). These findings suggest two rules for single crystal phases of high-T c hydrides: (I) phonons affect dT c /dp mainly through the parameter λ; and (II) dλ/dpis primarily determined by dζ/dp,ratherthandβ/dp. In view of ζbeing proportional to NϵF  ,thepressuretrendof NϵF  in LaH 9.63 (Fig. 4a) suggests a sign change of dT c /dpat studied pressures. Since the quantum fluxional nature of LaH 9.63 forbids a direct calculation of the Eliashberg function αωðÞ 2FωðÞby density functional methods, a quantitative confirmation of this is unachievable yet. However, considering that the change is so significant that even qualitative estimations could provide insights, we study the sign of dT c /dpslope in LaH 9.63 under two approximations, as the first step to access the problem: (I) calculating ωlog and  ω2based on the phonon spectrum FωðÞof LaH 9.63 , combined with αωðÞ 2approximated by that of quantum fcc-LaH 10 ; (II) approximating the βin LaH 9.63 by that of the quantum fcc-LaH 10 . The rationale of such approximation is based on empirical rules, and the fact that quantum fluxional LaH 9.63 largely retains the same local atomic environment as quantum fcc-LaH 10 (Supplementary Fig. 6). The FωðÞderived from the Fourier transform of the velocity autocorrelation functions in the CMD simulations at 240 K, as well as the approximated αωðÞ 2and αωðÞ 2FωðÞare shown in Supplementary Fig. 7. Withthisapproach,weevaluatedthepressuretrendofT c in LaH 9.63 using the AD formula (with parameters listed in Supplementary Table 2). It does indeed exhibit a positive dT c /dpslope at 137–163 GPa, in qualitative agreement with experimental observation1(Fig. 4b). Moreover, a negative dT c /dpslope is obtained at higher pressures in line with both experimental observations and theoretical results for fccLaH 10 (Supplementary Fig. 1). The standard deviation of NϵF  (Supplementary Fig. 5) reveals fluctuations of the electron structure near the Fermi level, which may affect the T c through parameters ζand λ. However, considering that superconductivity of LaH 9.6 is measured at a time scale far beyond that of the simulation, it is feasible to calculate the T c by a averaged NϵF  . Recently, the C2/mphase has been successfully prepared at 120 GPa by abrupt decompression of the fcc phase5.A positive dT c /dpslope measured in the subsequent compression was attributed to a boost in the T c due to phonon softening. As illustrated in the literature12,16, the softening of phonons results in a negative dT c /dp slope for fcc-LaH 10 , in contrast to the experimental observations below 150 GPa1. However, the concept of ‘ahigherT c near structural instability’5agrees well with our finding of a premature distortion of hydrogen substructure relative to an average fcc-La substructure in quantum fluxional LaH 9.63 , upon decompression to 137 GPa. In summary, the present work illustrates in lanthanum polyhydride the coexistence of quantum proton fluxion and hydrogen-induced highT c superconductivity. A premature distortion of hydrogen substructure upon volume expansion and its impact on the superconductivity through modulating the electronic structure is revealed. The findings 34.2 33.3 32.4 31.5 10-9 10-8 10-7 10-6 10-5 34.2 33.3 32.4 31.5 -0.0005 0.0000 0.0005 0.0010 0.0015 0.0020 0.0025 Hydrogen phase IV, 250 GPa Hydrogen phase IV, 350 GPa 240 K run1 240 K run2 120 K 60 K 240 K FeH, 202 GPa Cu2H, 96 GPa FeH, 33 GPa D (cm2/s) Volume ( ) Quantum: Classical: Cu2H, 15 GPa LaH9.63 137 143 150 163 176 Pressure (GPa) Pressure (GPa) XRD H | La |QF: QFS vs. FCC |QC: QFS vs. C2/m |QP: QFS vs. P1 (arb. units) Volume ( ) Decompression XRD 137 143 150 163 176 Å3Å3 ξ ba Fig. 3 | Pressure effects on the diffusivity and structural distortions. a The proton diffusion coefficient Dderived from the MSD in CMD simulations (run1) of 4 ps, and from the average MSD of four MD simulations of 24 ps. The Din three additional CMD simulations of 4 ps at 150GPa are also shown. These are two simulations at 120K and 60 K, and a simulation using higher total-energy convergence criteria (run2). The MSD in the CMD simulations is derived from the centroid trajectories. The measured Din Cu 2 H29 and FeH31, and calculated Din phase IV hydrogen30 are shown for comparison. bThe configurational distance ξ (with error bar indicating the standard deviation) of QFS from the crystallattices of staticfcc and C2mphases, and a triclinic P1 structure for H and La substructure. The P1 structure is built by scaling the lattice parameter a of a QFS sampled from the CMD simulation of LaH 9.63 at 176 GPa. Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 4 support the experimental argument concerning the hypostoichiometric nature of the superconducting samples, reveal the crucial role of vacancy defects in the structure, and provide a theoretical explanation for the positive dT c /dpslope of lanthanum polyhydride, with implications for understanding the structure and superconductivity of other hydrogen-rich materials. With significant progress of nuclear magnetic resonance techniques in diamond anvil cells, the measurement of proton mobility in metal hydrides is currently accessible at megabar pressures29,31.Weanticipatethatcontinuous experimental and theoretical studies of quantum fluxional structure of high-T c hydrides will stimulate new theoretical tools for understanding the superconducting behavior of quantum fluxional materials that cannot be precisely described by a single underlying static structure. Methods The ab initio calculations were performed using the plane-wave pseudopotential method, as implemented in the Vienna ab initio simulation program (VASP)36,withthePerdew–Burke–Ernzerhof (PBE) Generalized Gradient Approximation (GGA) density functional37,and the bare ion Coulomb potential treated in the projector augmented wave (PAW) framework38. The vacancy formation enthalpy of fcc-LaH 10 wascalculated from vacancy structures of fcc-LaH 10 and B 2 /n-H 2 39,with enthalpies calculated using the third-order Birch–Murnaghan isothermal equation of state (EOS)40 from ab initio energy-volume relations, asimplementedinthe EOS code41. The vacancymigrationenergy wascalculated by the climbing image nudged elastic band method (CINEB)42. The quantum nuclear dynamics were studied by path-integral molecular dynamics (PIMD) and centroid molecular dynamics (CMD)30, as implement in the PIMD code43. The superconductivity was analyzed using McMillan’stheory 34 with the T c estimated by the AllenDynes equations35. The crystal structures and MD trajectories were visualized by the VESTA44 and OVITO softwares45 respectively. More details are shown in the Supplementary Information. Data availability All data are available in the paper and from the author upon request. References 1. Drozdov, A. P. et al. Superconductivity at 250 K in lanthanum hydride under high pressures. Nature 569,528–531 (2019). 2. Somayazulu, M. et al. Evidence for superconductivity above 260 K in lanthanum superhydride at megabar pressures. Phys.Rev.Lett. 122, 027001 (2019). 3. Geballe, Z. M. et al. Synthesis and stability of lanthanum superhydrides. Angew. Chem. Int. Ed. 57,688–692 (2018). 4. Struzhkin, V. et al. Superconductivity in La and Y hydrides: remaining questions to experiment and theory. Mat. Rad. Extr. 5, 028201 (2020). 5. Sun, D. et al. High-temperature superconductivity on the verge of a structural instability in lanthanum superhydride. Nat. Commun. 12, 6863 (2021). 6. Snider, E. et al. Synthesis of yttrium superhydride superconductor with a transition temperature up to 262 K by catalytic hydrogenation at high pressures. Phys. Rev. Lett. 126, 117003 (2021). 7. Kong, P. et al. Superconductivity up to 243 K in the yttrium-hydrogen system under high pressure. Nat. Commun. 12, 5075 (2021). 8. Flores-Livas, J. A. et al. A perspective on conventional hightemperature superconductors at high pressure: Methods and materials. Phys. Rep. 856,1–78 (2020). 9. Pickard, C. J., Errea, I. & Eremets, M. I. Superconducting hydrides under pressure. Annu. Rev. Condens. Matter Phys. 11, 57–76 (2020). 10. Bardeen, J., Cooper, L. N. & Schrieffer, J. R. Theory of superconductivity. Phys. Rev. 108, 1175–1204 (1957). 11. Peng, F. et al. Hydrogen clathrate structures in rare earth hydrides at high pressures: possible route to room-temperature superconductivity. Phys. Rev. Lett. 119, 107001 (2017). Fig. 4 | Pressure effects on electronic properties and superconductivity. a The pressure trend of the averaged electronic density of states at the Fermi level NϵF  in quantum fluxional LaH 9.63 (with distributions and standard deviation of statics shown in Supplementary Fig. 5), and the pressure trend of NϵF  in static P1-LaH 9.63 , fcc-LaH 10 ,andC2/m-LaH 10 .TheNϵF  of LaH 9.63 in fcc and C2/mphase estimated using therigid-band approximation are also shown for comparison. bThe pressure trend of T c in quantum fluxional LaH 9.63 calculated based on various Gaussian smearings (σ) of the phonon spectrum FωðÞ, together with the values measured for LaH 9.6 by Drozdov et al.1and those derived from AD, Migdal–Eliashberg equations (ME), and superconducting DFT (SCDFT) for fcc-LaH 10 by Errea et al.16. Fitted lines are shown to guide the sight. Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 5 12. Liu, H., Naumov, I. I., Hoffmann, R., Ashcroft, N. W. & Hemley, R. J. Potential high-T c superconducting lanthanum and yttrium hydrides at high pressure. Proc.NatlAcad.Sci.USA114, 6990 (2017). 13. Liu, L. et al. Microscopic mechanism of room-temperature superconductivity in compressed LaH 10 .Phys.Rev.B99,140501(2019). 14. Wang, C., Yi, S. & Cho, J.-H. Pressure dependence of the superconducting transition temperature of compressed LaH 10 .Phys. Rev. B100,060502(2019). 15. Quan, Y., Ghosh, S. S. & Pickett, W. E. Compressed hydrides as metallic hydrogen superconductors. Phys. Rev. B 100, 184505 (2019). 16. Errea, I. et al. Quantum crystal structure in the 250-kelvin superconducting lanthanum hydride. Nature 578,66–69 (2020). 17. Kruglov, I. A. et al. Superconductivity of LaH 10 and LaH 16 polyhydrides. Phys.Rev.B101, 024508 (2020). 18. Wang,H.,Tse,J.S.,Tanaka,K.,Iitaka,T.&Ma,Y.Superconductive sodalite-like clathrate calcium hydride at high pressures. Proc. Natl Acad.Sci.USA109, 6463 (2012). 19. Xie, H. et al. High-temperature superconductivity in ternary clathrate YCaH 12 under high pressures. J. Phys: Condens Matter 31, 245404 (2019). 20. Song,H.etal.HighT c superconductivity in heavy rare earth hydrides. Chin. Phys. Lett. 38, 107401 (2021). 21. Ma, L. et al. High-temperature superconducting phase in clathrate calcium hydride CaH 6 up to 215 K at a pressure of 172 GPa. Phys. Rev. Lett. 128, 167001 (2022). 22. Goncharov, A. F., Chuvashova, I., Ji, C. & Mao, H.-k Intermolecular coupling and fluxional behavior of hydrogen in phase IV. Proc. Natl Acad.Sci.USA116, 25512 (2019). 23. Howie, R. T., Guillaume, C. L., Scheler, T., Goncharov, A. F. & Gregoryanz, E. Mixed molecular and atomic phase of dense hydrogen. Phys.Rev.Lett.108,125501(2012). 24. Zha, C. et al. High-pressure measurements of hydrogen phase IV using synchrotron infrared spectroscopy. Phys. Rev. Lett. 110, 217402 (2013). 25. Magdău,I.B.&Ackland,G.J.Identification of high-pressure phases III and IV in hydrogen: simulating Raman spectra using molecular dynamics. Phys. Rev. B 87, 174110 (2013). 26. Shiga, M., Tachikawa, M. & Miura, S. A unified scheme for ab initio molecular orbital theory and path integral molecular dynamics. J. Chem. Phys. 115,9149–9159 (2001). 27. Liu, H. et al. Dynamics and superconductivity in compressed lanthanum superhydride. Phys. Rev. B 98, 100102(R) (2018). 28. Katoh,E.,Yamawaki,H.,Fujihisa,H.,Sakashita,M.&Aoki,K.Protonic diffusion in high-pressure ice VII. Science 295, 1264 (2002). 29. Meier, T. et al. Proton mobility in metallic copper hydride from highpressure nuclear magnetic resonance. Phys.Rev.B102,165109 (2020). 30. Liu,H.&Ma,Y.ProtonordeuterontransferinphaseIVofsolid hydrogen and deuterium. Phys. Rev. Lett. 110, 025903 (2013). 31. Meier, T. et al. Pressure-induced hydrogen-hydrogen interaction in metallic FeH revealed by NMR. Phys.Rev.X9,031008(2019). 32. Zhu, L. et al. A fingerprint based metric for measuring similarities of crystalline structures. J. Chem. Phys. 144, 034203 (2016). 33. Boeri L. in Handbook of Materials Modeling: Applications: Current and Emerging Materials (eds Andreoni, W. & Yip, S.) (Springer International Publishing, 2018). 34. McMillan, W. L. Transition temperature of strong-coupled superconductors. Phys. Rev. 167, 331–344 (1968). 35. Allen, P. B. & Dynes, R. C. Transition temperature of strong-coupled superconductors reanalyzed. Phys. Rev. B 12,905–922 (1975). 36. Kresse, G. & Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys.Rev.B 54, 11169–11186 (1996). 37. Perdew, J. P., Burke, K. & Ernzerhof, M. Generalized gradient approximation made simple. Phys.Rev.Lett.77, 3865–3868 (1996). 38. Blöchl, P. E. Projector augmented-wave method. Phys.Rev.B50, 17953–17979 (1994). 39. Pickard, C. J. & Needs, R. J. Structure of phase III of solid hydrogen. Nat. Phys. 3,473–476 (2007). 40. Birch, F. Finite elastic strain of cubic crystals. Phys. Rev. 71, 809–824 (1947). 41. Dewhurst, J. & Ambrosch-Draxl C. The EXCITING Code Users’ Manual Version 0.9.74 (2021). 42. Henkelman,G.,Uberuaga,B.P.&Jónsson,H.Aclimbingimage nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 113,9901–9904 (2000). 43. Shiga, M., Tachikawa, M. & Miura, S. Ab initio molecular orbital calculation considering the quantum mechanical effect of nuclei by path integral molecular dynamics. Chem. Phys. Lett. 332, 396–402 (2000). 44. Momma, K. & Izumi, F. VESTA3 for three-dimensional visualization of crystal, volumetric and morphology data. J. Appl. Crystallogr. 44, 1272–1276 (2011). 45. Stukowski, A. Visualization and analysis of atomistic simulation data with OVITO-the Open Visualization Tool. Model. Simul. Mater. Sci. Eng. 18, 015012 (2010). Acknowledgements H.W. is thankful to Y. Ma and M. Shiga for valuable discussions, to V. Minkov and M. Eremets for the XRD data, and to C.Wang for the EPC data. The project is supported by the National Natural Science Foundation of China (Grant No. 11974135, 11874176, 12174170, and 12074138), the Natural Sciences and Engineering Research Council of Canada, the EPSRC through grants EP/P022596/1, and EP/S021981/1, and the startup funds of the office of the Dean of SASN of Rutgers University-Newark. P. T. S. thanks the Department of Materials Science and Metallurgy at the University of Cambridge for generous funding. The work of P. T. S. is further supported through a Trinity Hall research studentship. I. E. acknowledges financial support by the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program (grant agreement no. 802533). We used the computing facilities at Beijing Super Cloud Computing Center. Author contributions H.W., P.S., I.E., F.P., Z.L., H.L., and L.Z. performed the calculations. H.W., Y.Y., P.S., and C.P. wrote the manuscript, with input from all co-authors. Competing interests The authors declare no competing interests. Additional information Supplementary information The online version contains supplementary material available at https://doi.org/10.1038/s41467-023-37295-1. Correspondence and requests for materials should be addressed to Hui Wang. Peer review information Nature Communications thanks the anonymous reviewers for their contribution to the peer review of this work. Peer reviewer reports are available. Reprints and permissions information is available at http://www.nature.com/reprints Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations. Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 6 Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/ licenses/by/4.0/. © The Author(s) 2023 Article https://doi.org/10.1038/s41467-023-37295-1 Nature Communications | (2023) 14:1674 7 Supplementary Information of “Quantum structural fluxion in superconducting lanthanum polyhydride” Hui Wang1,2*, Pascal T. Salzbrenner3, Ion Errea4,5,6, Feng Peng7, Ziheng Lu3, Hanyu Liu2,8, Li Zhu9, Chris J. Pickard3,10 and Yansun Yao11 1Key Laboratory for Photonic and Electronic Bandgap Materials (Ministry of Education), School of Physics and Electronic Engineering, Harbin Normal University, Harbin 150025, China 2International Center for Computational Method & Software, College of Physics, Jilin University, Changchun 130012, China 3Department of Materials Science & Metallurgy, University of Cambridge, 27 Charles Babbage Road, Cambridge CB3 0FS, United Kingdom 4Fisika Aplikatua Saila, Gipuzkoako Ingeniaritza Eskola, University of the Basque Country (UPV/EHU), Europa Plaza 1, 20018 Donostia/San Sebastián, Spain 5Centro de Física de Materiales (CSIC-UPV/EHU), Manuel de Lardizabal Pasealekua 5, 20018 Donostia/San Sebastián, Spain 6Donostia International Physics Center (DIPC), Manuel de Lardizabal Pasealekua 4, 20018 Donostia/San Sebastián, Spain 7College of Physics and Electronic Information, Luoyang Normal University, Luoyang 471022, P. R. China 8State Key Laboratory of Superhard Materials and International Center of Future Science, Jilin University, Changchun 130012, China 9Department of Physics, Rutgers University, Newark, NJ 07102, USA 10Advanced Institute for Materials Research, Tohoku University 2-1-1 Katahira, Aoba, Sendai, 980-8577, Japan 11Department of Physics and Engineering Physics, University of Saskatchewan, Saskatoon, Saskatchewan S7N 5E2, Canada *e-mail: [email protected] Table of contents 1. Extended Data Figures; Page 1 ~ 7 Supplementary Figure 1 ~ Supplementary Figure 7 2. Extended Data Tables; Page 8 ~ 9 Supplementary Table 1 ~ Supplementary Table 2 3. Computational Details; Page 10 ~ 15 4. Supplementary References; Page 16 ~ 17 Supplementary Figure 1 Supplementary Figure 1| The ‘positive dTc/dp contradiction’. Recent experimental studies by Drozdov et. al.1, Somayazulu et. al.2, Struzhkin et. al.3, Kong et. al.4, Snider et. al.5 and Ma et. al.6 discovered several near-room temperature superconducting polyhydrides, and the measured pressure trend of Tc usually reveals a positive dTc/dp at some pressure regimes, in contradiction with the negative slope predicted by earlier or subsequent calculations by Peng et. al.7, Liu H. et. al.8, Liu L. et. al.9, Wang C. et. al.10, Quan et. al.11, Errea et. al.12, Kruglov et. al.13, Wang H. et. al.14 and Duan et. al.15, 16 on the same hydride based on BCS theory. The contradiction is clearly seen from the data. The present work focuses on the 250-kelvin lanthanum polyhydride, which has a positive dTc/dp measured by Drozdov et. al. on LaH9.61, and a negative dTc/dp calculated by Errea et. al.12 on fcc-LaH10 at pressures of 137-150 GPa. 1／17 Supplementary Table 1 Supplementary Table 1| Parameters impacting the pressure trend of Tc. Based on retrievable Ω(p) (where Ω represents 𝜔𝑙𝑜𝑔 or 𝜔2) (Ref 11) or 𝛼(𝜔)2𝐹(𝜔) from the literature (Ref 12, 8, 10 and 14), we derived dΩ/dp and analyzed the role played by several parameters in the pressure trend of Tc using AD equations, with the Coulomb coupling constant μ* set to 0.1. The resulting pressure trends of Tc for cubic LaH10 and CaH6 are in good agreement with the literature (Supplementary Figure 1). The analysis demonstrates that 𝜁 decreases more than two times faster than 𝛽 increases upon compression, which leads to a negative sign of the dλ/dp slope. This suggests that d𝜁/dp dominates the pressure trend of λ. λ has a monotonously negative pressure dependence, but decreases less significantly than 𝜔𝑙𝑜𝑔 increases. However, due to the nonlinear dependence of Tc on λ, either through the well-known exponential function or the ‘strong-coupling correction’ and ‘shape correction’ factors (see details in page 14), dλ/dp plays a dominant role in determining the dTc/dp slope. Our results for LaH9.63 are also shown in order to enable a direct comparison. Our d𝜔𝑙𝑜𝑔 /dp is comparable to those of cubic LaH10 (fcc) and CaH6 (bcc), but d𝜁/dp decreases significantly at lower pressures, which reverses the sign of the dTc/dp slope through a considerable reduction of dλ/dp below 163 GPa. We further estimated the limit of stability of the +dTc/dp slope against variations of 𝛽 and d 𝛽 /dp: I) at a fixed d 𝛽 /dp of 0.01235 eV3/megabar, decreasing (increasing) 𝛽 by 50 % results in a decrease (an increase) of dTc/dp by 59% (41%), but does not change the positive sign; II) at a fixed 𝛽 of 0.04352 eV3 at 176 GPa, an about 127 % decrease of the d 𝛽 /dp slope is needed to invert the sign of dTc/dp. Literature data generally report a positive d 𝛽 /dp slope for cubic hydrides covering a range of H/M ratios (M= La and Ca) from 10 to 6. Therefore, LaH9.63 is very unlikely to exhibit a -d 𝛽 /dp slop in qualitative contradiction with most other hydrides. The analysis indicates that the +dTc/dp slope is robust to errors induced by the approximation of 𝛽 and d 𝛽 /dp in LaH9.63 by those in quantum fcc-LaH10. 8／17 Supplementary Table 2 Supplementary Table 2| Parameters for estimating the pressure trend of Tc. Parameters for calculating the pressure trend of Tc shown in Fig 4b with smearing 𝜎 = 0.0025 eV for F(ω) are listed here for illustration. Furthermore, we compare the pressure trends of Tc calculated for various 𝜎 of 𝛼(𝜔)2 (i.e. 0, 0.01 and 0.05). It is found that the trend is determined by the overall shape of 𝛼(𝜔)2. Due to the inexistence of a method that can directly calculate the 𝛼2𝐹(𝜔) of a quantum fluxional structure, we estimated frequency moments Ω of 𝛼2𝐹(𝜔) in LaH9.63 based on the calculated 𝐹(𝜔) of LaH9.63 and the 𝛼(𝜔)2 approximated by that of quantum fcc-LaH1012. To evaluate the potential errors associated with the approximation, we calculated the Tc using reshaped 𝛼2𝐹(𝜔) obtained by artificial scaling of 𝛥𝜔. Empirically, EPC weights the phonon spectrum to lower frequencies in high-Tc hydrides. This weighting effect changes smoothly with pressure, as observed from the shape of 𝛼(𝜔)2 in quantum fcc-LaH10 (Supplementary Figure 7). Considering this, and that quantum LaH9.63 and fcc-LaH10 exhibit a similar phonon hardening trend with increasing pressure as well as similar shapes of 𝐹(𝜔), our results suggest that the pressure trend of Tc is robust to variations of Ω. This indicates that it is practicable to estimate qualitatively the sign of dTc/dp slop in LaH9.63 under this approximation. 9／17 Computational Details Molecular dynamics calculations Quantum nuclear dynamics were studied using path-integral molecular dynamics (PIMD) and centroid molecular dynamics (CMD), with the massive Nosé-Hoover chain (NHC) thermostats on NVT and NVE ensembles (N-number of particles; V-volume; T-temperature, and E-energy), respectively, as implemented in the PIMD code17. NpT-PIMD (ppressure) simulations were performed to estimate the quantum pressure-volume relations. Beads number was set to 16 for the cases without a specification, and time step was set to 0.5 and 0.05 fs for PIMD and CMD, respectively. A short PIMD simulation was always carried out on the initial structure to prepare the pre-equilibrium state in order to promote the efficiency and stability of the CMD simulation. The parameters of the standard ab initio MD simulations are the same as those of the PIMD simulations, except for a reduced bead number of 1. The MD runs of LaH9.63 were initialized from different centroid configurations of the 4-picosecond CMD trajectory with a sampling interval of 0.25 ps for 16 short runs of 12000 steps and of 1 ps for 4 long runs of 52000 steps. The underlying ab initio total energy and forces were calculated based on the plane-wave pseudopotential method, as implemented in the Vienna ab initio simulation program (VASP)18. The convergence criterion for the total energy (Ecc) was chosen to be 3×10−6 eV/atom with a plane-wave cutoff of 325 eV and Γ-point sampling of the first Brillouin zone. Decreasing Ecc to 3×10−7 eV/atom (run 2) results in a nearly negligible change of the proton diffusion coefficient D in LaH9.63 at 150 GPa (see Fig. 3a in the main text). The volume of the cubic simulation cells of LaH9.63 were fixed at 1097.70, 1079.98, 1062.39, 1038.16, and 1014.35 Å3. The resulting temperature fluctuations at equilibrium are within ± 15 K, ± 10 K and ± 25 K with respect to the target in the PIMD, CMD and MD simulations, respectively. Other details are summarized in Supplementary Table 3. Supplementary Table 3. Details of simulations. Phase Number of Atoms, H : La Pressures, GPa Temperature, K Method Steps runs LaH9.63 308:32 150 60, 120, 240 MD 12000 1-16 137, 143, 150, 163,176 240 MD 52000 1-4 150 240 PIMD 2000 1,2 CMD 82000 1,2 137, 143, 163,176 240 PIMD 2000 1 CMD 82000 1 150 60 PIMD 8000 3§ 60, 120 PIMD 2000 1 CMD 82000 1 LaH9.0 288:32 130, 137, 150, 170, 180 300 PIMD† 2000 1 LaH10 320:32 130, 137, 150, 170, 180 300 PIMD† 2000 1 137, 150 240 PIMD 2000 1 CMD 42000 1 § NVT simulations with 16, 32 and 64 beads were performed to check the impact of bead number on MSD; †NpT simulations. 10／17 The quantum pressure-volume relation of fcc-LaH10 and fcc-LaH9. The quantum pressure-volume relation of fcc-LaH10 and fcc-LaH9 at 300 K was estimated from the volume averaged over 1 ps at equilibrium obtained in NpT-PIMD simulations using a cubic box and Martyna’s equation of motion, as implemented in the PIMD code17. Our results reveal a H/La ratio of 9.54 (or 9.71) for the experimental samples synthesized at 150 GPa (or 137 GPa), in good agreement with experimental estimation of a LaH9.61. To ensure our pressure trend of Tc can directly compare with the experimental one, we built structural models at 137 and 150 GPa using experimental lattice parameters. This guarantees that the coincidence of this pressure range (i.e. ‘137-150 GPa’) between calculation and experiments does not dependent on the experimental or theoretical pressure scales. Notably, this is indeed the pressure range that the positive dTc/dp slope were measured experimentally and calculated in this work. The PBE functionals may lead to potential shift of the calculated Tc v.s. pressure curve in relative to the experimental measured one at pressure above 150 GPa, however, which does not affect our conclusion. The quantum pressure of 143, 163 and 176 GPa were estimated based on interpolation. The mean square displacement, 〈∆𝒓𝟐〉 〈∆𝑟(𝑡)𝟐〉 = 〈|𝑟𝑖(𝑡)−𝑟0(𝑡)−[𝑟𝑐𝑚(𝑡)−𝑟𝑐𝑚(0)]|2〉𝐼, where 𝑟𝑖(𝑡) is the position of the diffusive proton i, which exhibits at least one jump over a threshold value of 0.7 Å over the course of the simulation, 𝑟𝑐𝑚(𝑡) represents the position of the center of mass of the system at the time t, and the 〈 〉 represent an average over 10 time steps and the total number of ‘particle’ I, which equals 308 for proton and 12 for vacancy in LaH9.63. In the simulations, some of diffusive protons hop back to the previously occupied position occasionally due to ‘traffic jam’, which results in a decrease of 〈∆𝒓𝟐〉 . The proton traffic jam and local structural relaxation result in plateaus in 〈∆𝒓𝟐〉 curves in both the CMD and MD simulations. However, 〈∆𝒓𝟐〉 generally increases with increasing of t, which is clearly revealed by long-time MD simulations of 24 ps. Supplementary Figure 8. MSD of PIMD simulation with various bead numbers at 60 K. MSDs derived from the centroid trajectory of PIMD simulations at 150 GPa and an equilibrium temperature of 60 K with 16, 32 and 64 beads are shown in Supplementary Figure 8. All curves suggest quantum proton diffusion at this temperature. Considering that the strength of NQE increases with decreasing temperature, larger bead number gives better MSD. The difference of 〈∆𝒓𝟐〉 values at the end of the simulation (i.e. 2 ps) is not significant, which may be attributed to the ‘traffic jam’ of protons in LaH9.63 associated with the small amount of stoichiometric defect. 11／17 The diffusion coefficient, D We approximate the proton diffusion coefficient D in LaH9.63 from the slope of 〈∆𝑟𝟐〉 obtained from a linear fit, as shown in Supplementary Figure 9 for CMD simulations run1 and run2 at 150 GPa and 240 K. D is computed from 𝐷 = 𝑠𝑙𝑜𝑝𝑒(〈∆𝑟𝟐〉) 6 ⁄. Supplementary Figure 9. Linear fit of MSD. Vacancy formation enthalpy, Hf We carried out a full variable-cell optimization (using ISIF=3 in VASP) of fcc-LaH107, 8 (in a 2×2×2 supercell of the conventional unit cell containing 352 atoms) and B2/n-H2 (in the conventional unit cell of 48 atoms)19 at pressures between 120 and 220 GPa with an interval of 5 GPa. We employed k-mesh grids of sizes 3×3×3 and 15×9×9, respectively, to carry out the required DFT calculations. The VT and VC LaH9.97 vacancy structures were then constructed by removing the appropriate H atom from the relaxed fcc-LaH10 structure at each pressure, as described in the main text. Variable-cell optimization of these structures induced proton diffusion, resulting in a distortion of the structures. In order to isolate the effect of the vacancies, therefore, we only relaxed the atomic positions, keeping the cell shape and volume fixed (ISIF=2 in VASP). We used a k-mesh grid equivalent to that for fcc-LaH10. A cut-off energy of 325 eV and a convergence criterion for the total energy of 3×10−7 eV/atom were employed in all calculations. Using the resulting V-E relationship of the relaxed structures, the pressure-enthalpy curves were then calculated using the third-order Birch–Murnaghan isothermal equation of state20. From the enthalpies for each pressure, the Hf of the VC and VT structures were calculated as follows, Hf = H(LaH9.97) – H(fcc-LaH10) + H(B2/n-H2)/48, with the configurationally averaged Hf expressed as, Hf(average) = 0.2Hf(VT) + 0.8Hf(VC). The calculated pressure-volume curve gives a pressure shift of ~5 GPa for LaH9.97 with respect to the experimental measurements on fcc-LaH9.6 at 120-220 GPa, which does not affect our discussions. It should be noted that Hf is calculated classically in the static lattice approximation on the Born–Oppenheimer energy surface, which differs from the configurational energy surface of the quantum crystal (e.g., fcc-LaH10 in Ref. 12) or the free energy surface of the classical or quantum crystal at finite temperature. The stability and pressure boundary of the 12／17 vacancy structure may be affected by more complex physical effects existing in the superconducting samples, for instance nuclear quantum effects (including zero-point motion of nuclei and proton tunnelling), or temperature effects. Indeed, temperature effects only enhance the formation of vacancies. In this sense, the conclusion from the Hf curves shown in Fig.1a of the main text, namely that vacancy formation is favored below 158 GPa in the classical crystal picture, while being based on a model calculation subject to many approximations, indicates that lanthanum superhydride should be rather prone to vacancy formation. The configurational distance, ξ The crystal fingerprint technique utilized in this study is based on Gaussian overlap matrices, which represent the local environments of all atoms in a unit cell and can efficiently determine configurational distances ξ between crystalline structures satisfying the mathematical requirements of a metric21. The ξ of the QFS of LaH9.63 vs. various static structures (e.g. fcc, C2/m and P1) were obtained by statistical analysis over the ξ of 400 centroid configurations sampled from 16bead CMD trajectories of 4 ps with a time interval of 10 fs at each pressure. Electronic density of states at the Fermi level, 𝑵(𝝐𝑭) The 𝑁(𝜖𝐹) is of central importance to the understanding of conventional superconductivity for hydrogen-rich materials. It is quite general that the pressure dependence of 𝑁(𝜖𝐹) directly correlates to the pressure trend of the superconducting Tc. The 𝑁(𝜖𝐹) of the quantum structures of LaH9.63 was calculated by carrying out a statistical average over the 𝑁(𝜖𝐹) of 400 centroid configurations sampled from 16-bead CMD trajectories of 4 ps with a time interval of 10 fs at each pressure. The 𝑁(𝜖𝐹) of each centroid configuration was calculated using the tetrahedron method with Blöchl corrections22, with a convergence criterion for the total energy of 3×10−6 eV/atom for a cut-off energy of 325 eV and a 5×5×5 k-mesh grid. The 𝑁(𝜖𝐹) of static fcc-LaH10, C2/m-LaH10, and P1-LaH10 were calculated based on primitive cells by the tetrahedron method with Blöchl corrections, with a convergence criterion for the total energy of 3×10−6 eV/atom for a cut-off energy of 325 eV and 21×21×21, 27×27×15 and 5×5×5 k-mesh grids, respectively. For the fcc and C2/m phases, the 𝑁(𝜖𝐹) of static LaH9.63 were estimated based on the rigid-band model23 of the electronic structure of LaH10 by virtue of an artificial shift of 𝜖𝐹. Supplementary Figure 10. The 𝑁(𝜖𝐹) and 𝜉 (of H substructure) in LaH9.63 compared to those of LaH10. The f.u. is the abbreviation of formula unit, throughout the article. 13／17 Comparison on 𝑵(𝝐𝑭) and ξ between LaH9.6 and LaH10 As shown in Supplementary Figure 10, we compared the configurationally averaged 𝜉 and 𝑁(𝜖𝐹) in quantum LaH9.63 and LaH10 in the pressure range of 137–150 GPa. The results suggest that the H sublattice distortion is suppressed in the vacancy-free case. Moreover, a negative d𝑁(𝜖𝐹)/dp slop is observed in LaH10, in contrast to the positive d𝑁(𝜖𝐹)/dp slop of LaH9.63. The frequency moments of 𝜶𝟐𝑭(𝝎) in LaH9.63, 𝝎𝒍𝒐𝒈 and 𝝎 𝟐 The 𝛼2𝐹(𝜔) of LaH9.63 is estimated based on the 𝐹(𝜔) of LaH9.63 derived by Fourier transforming the velocity autocorrelation functions (VACF) in the CMD simulations at 240 K in combination with coupling functions 𝛼(𝜔)2 approximated by that of quantum fcc-LaH10: 1. Deriving the 𝐹(𝜔) of LaH9.63 by Fourier transforming the VACF in the CMD simulations at 240 K and volume (V) of 34.30, 33.75, 33.20, 32.44, and 31.70 Å3/f.u.; 2. Calculating the 𝛼(𝜔)2 from 𝛼2𝐹(𝜔) and 𝐹(𝜔) in quantum fcc-LaH10: 𝛼(𝜔)2= 𝛼2𝐹(𝜔)/ 𝐹(𝜔) at volume (V′) of 35.23, 32.97 and 30.36 Å3/f.u., and deriving the 𝛼(𝜔)2 at V by inverse volume weight interpolation of 𝛼(𝜔)2 of neighboring V′; 3. Calculating the 𝛼(𝜔)2𝐹(𝜔) of LaH9.63 based on 𝐹(𝜔) of LaH9.63 and 𝛼(𝜔)2 of quantum fcc-LaH10: 𝛼2𝐹(𝜔)_LaH9.63 = 𝛼(𝜔)2_LaH10 × 𝐹(𝜔)_LaH9.63 at various V values; 4. Calculating the 𝜔𝑙𝑜𝑔 and 𝜔2 following their definitions in Ref. 24. The frequency step 𝛥𝜔 is set to 0.0004 eV in the calculations. To include the vibrational modes in both of lowand high-frequency regime, the 𝐹(𝜔) are derived from the VACFs of the La and H sublattices separately. Since smeared 𝐹(𝜔) is generally adopted to converge EPC parameters, we smear 𝐹(𝜔) with Gaussian profiles with widths (𝜎) of 0.002 ~ 0.003 eV. 𝜔𝑙𝑜𝑔 and 𝜔2 are therefore obtained following their definitions in Ref. 24. In addition, we reshaped the 𝛼2𝐹(𝜔) of LaH9.63 by artificial scaling of 𝛥𝜔 to 0.0003 and 0.0005 eV to evaluate the potential errors associated with the approximation in step 3, see the caption of Supplementary Table 2. The electron-phonon coupling constant, λ We estimate λ in LaH9.63 using McMillan’s formula25 (𝜆 = 𝑁(𝜖𝐹)⟨𝛪2⟩/𝑀𝜔2 2) as follows: 1. Combining variables as 𝜆 = 𝛽 ⋅𝜁 with 𝛽 = ⟨𝛪2⟩/𝑀 and 𝜁 = 𝑁(𝜖𝐹)/𝜔2 2; 2. Estimating 𝜁(𝑉) based on 𝑁(𝜖𝐹) and 𝜔2 of LaH9.63 at V of 34.30, 33.75, 33.20, 32.44, and 31.70 Å3/f.u.; 3. Calculating 𝛽(V′) from 𝜆, 𝑁(𝜖𝐹) and 𝜔2of quantum fcc-LaH10 in Ref. 12 at V′ of 35.23, 32.97 and 30.36 Å3/f.u. and deriving the 𝛽(V) by inverse volume weight interpolation of neighboring 𝛽(V′); 4. Calculating 𝜆(𝑉) in LaH9.63 using the estimated 𝜁(𝑉) of LaH9.63 and the 𝛽(V) of quantum fcc-LaH10. Critical temperature of superconductivity, Tc The Tc was evaluated using the Allen–Dynes-modified McMillan equation (AD)24, 𝑇𝑐=𝑓1𝑓2𝜔𝑙𝑜𝑔 1.2 𝑒𝑥𝑝[−1.04(1+𝜆) 𝜆−𝜇∗(1+0.62𝜆)], where 𝑓1 and 𝑓2 are the ‘strong-coupling correction’ and ‘shape correction’ factors for strong-coupling systems: 𝑓1= [1+(𝜆/  1)3/2]1/3 with  1= 2.46(1+3.8  ∗); 𝑓2= 1+(𝜔2𝜔𝑙𝑜𝑔 ⁄−1)𝜆2(𝜆2+  2 2) ⁄ with  2= 1.82(1+6.3  ∗)(𝜔2𝜔𝑙𝑜𝑔 ⁄). μ* is the Coulomb coupling constant, which is set to a typical value of 0.1 in this wok. 14／17 Structural information of static LaH10-𝛿. The simulation cells of fcc-LaH9 (320 atoms) and fcc-LaH10 (352 atoms) are initialized by a 2×2×2 extension of the unit cell of F-43m-LaH913 and Fm-3m-LaH1012, respectively. We initialize the LaH9.63 structure (340 atoms) by adding randomly 20 H atoms to the simulation cell of fcc-LaH9, constraining neighboring vacancies to be no less than 4.4 Å apart to ensure an approximately uniform vacancy distribution, as shown in Supplementary Figure 11. A PIMD simulation (2000 steps) was always carried out on the initial structure to prepare the pre-equilibrium state for the CMD simulation. Vacancy diffusion generally results in a random distribution of vacancies on the 8c and 32f Wyckoff position in Fm-3m-LaH10 within several hundred PIMD steps. The static P1-LaH9.63, mimicking the quantum fluxional LaH9.63, is built by scaling the lattice parameter a of a QFS sampled from the CMD simulation of LaH9.63 at 176 GPa. The initial vacancy structure of LaH9.63 and static P1-LaH9.63 at 150 GPa can be found at the link— https://doi.org/10.6084/m9.figshare.20484420.v1. Supplementary Figure 11. Initial vacancy structure of LaH9.63. 15／17 Supplementary References (Repeating citations for some references of the main text) 1. Drozdov AP , et al. Superconductivity at 250 K in lanthanum hydride under high pressures. Nature 569, 528-531 (2019). 2. Somayazulu M , et al. Evidence for Superconductivity above 260 K in Lanthanum Superhydride at Megabar Pressures. Phys Rev Lett 122, 027001 (2019). 3. Struzhkin V , et al. Superconductivity in La and Y hydrides: Remaining questions to experiment and theory. Mat Rad Extr 5, 028201 (2020). 4. Kong P , et al. Superconductivity up to 243 K in the yttrium-hydrogen system under high pressure. Nat Commun 12, 5075 (2021). 5. Snider E , et al. Synthesis of Yttrium Superhydride Superconductor with a Transition Temperature up to 262 K by Catalytic Hydrogenation at High Pressures. Phys Rev Lett 126, 117003 (2021). 6. Ma L , et al. High-Temperature Superconducting Phase in Clathrate Calcium Hydride CaH6 up to 215 K at a Pressure of 172 GPa. Phys Rev Lett 128, 167001 (2022). 7. Peng F, Sun Y, Pickard CJ, Needs RJ, Wu Q, Ma Y. Hydrogen Clathrate Structures in Rare Earth Hydrides at High Pressures: Possible Route to Room-Temperature Superconductivity. Phys Rev Lett 119, 107001 (2017). 8. Liu H, Naumov II, Hoffmann R, Ashcroft NW, Hemley RJ. Potential high-Tc superconducting lanthanum and yttrium hydrides at high pressure. Proc Natl Acad Sci USA 114, 6990 (2017). 9. Liu L, Wang C, Yi S, Kim KW, Kim J, Cho J-H. Microscopic mechanism of room-temperature superconductivity in compressed LaH10. Phys Rev B 99, 140501 (2019). 10. Wang C, Yi S, Cho J-H. Pressure dependence of the superconducting transition temperature of compressed LaH10. Phys Rev B 100, 060502 (2019). 11. Quan Y, Ghosh SS, Pickett WE. Compressed hydrides as metallic hydrogen superconductors. Phys Rev B 100, 184505 (2019). 12. Errea I , et al. Quantum crystal structure in the 250-kelvin superconducting lanthanum hydride. Nature 578, 66-69 (2020). 13. Kruglov IA , et al. Superconductivity of LaH10 and LaH16 polyhydrides. Phys Rev B 101, 024508 (2020). 14. Wang H, Tse JS, Tanaka K, Iitaka T, Ma Y. Superconductive sodalite-like clathrate calcium hydride at high pressures. Proc Natl Acad Sci USA 109, 6463 (2012). 15. Xie H , et al. High-temperature superconductivity in ternary clathrate YCaH12 under high pressures. J Phys: Condens Matter 31, 245404 (2019). 16. Song H, Zhang Z, Cui T, Pickard CJ, Kresin VZ, Duan D. High Tc superconductivity in heavy Rare Earth Hydrides. Chin Phys Lett 38, 107401 (2021). 16／17 17. Shiga M, Tachikawa M, Miura S. A unified scheme for ab initio molecular orbital theory and path integral molecular dynamics. J Chem Phys 115, 9149-9159 (2001). 18. Kresse G, Furthmüller J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys Rev B 54, 11169-11186 (1996). 19. Pickard CJ, Needs RJ. Structure of phase III of solid hydrogen. Nat Phys 3, 473-476 (2007). 20. Birch F. Finite Elastic Strain of Cubic Crystals. Phys Rev 71, 809-824 (1947). 21. Zhu L , et al. A fingerprint based metric for measuring similarities of crystalline structures. J Chem Phys 144, 034203 (2016). 22. Blöchl PE. Projector augmented-wave method. Phys Rev B 50, 17953-17979 (1994). 23. Boeri L. Understanding Novel Superconductors with Ab Initio Calculations. In: Handbook of Materials Modeling: Applications: Current and Emerging Materials (eds Andreoni W, Yip S). Springer International Publishing (2018). 24. Allen PB, Dynes RC. Transition temperature of strong-coupled superconductors reanalyzed. Phys Rev B 12, 905922 (1975). 25. McMillan WL. Transition Temperature of Strong-Coupled Superconductors. Phys Rev 167, 331-344 (1968). 17／17