Excitonic insulator to superconductor phase transition in ultra-compressed helium

Helium, the second most abundant element in the universe, exhibits an extremely large electronic band gap of about 20 eV at ambient pressures. While the metallization pressure of helium has been accurately determined, thus far little attention has been paid to the specific mechanisms driving the band-gap closure and electronic properties of this quantum crystal in the terapascal regime (1 TPa = 10 Mbar). Here, we employ density functional theory and many-body perturbation calculations to fill up this knowledge gap. It is found that prior to reaching metallicity helium becomes an excitonic insulator (EI), an exotic state of matter in which electrostatically bound electron-hole pairs may form spontaneously. Furthermore, we predict metallic helium to be a superconductor with a critical temperature of ≈ 20 K just above its metallization pressure and of ≈ 70 K at 100 TPa. These unforeseen phenomena may be critical for improving our fundamental understanding and modeling of celestial bodies.

upon its closure is characteristic of a semimetal (i.e., the density of electronic states at the Fermi level is very small) 8 . Meanwhile, the lattice phonons involving atomic displacements perpendicular to the hcp basal plane drive the widening of the band gap at very high pressures 8 .
A detailed understanding of the electronic band structure properties of this archetypal quantum crystal 11 , however, is still lacking. First, about half a century ago the existence of an exotic insulating phase called "excitonic insulator" (EI) was predicted in which electrons and holes spontaneously form bound pairs called excitons 12 . The EI phase could be stabilized at sufficiently low temperatures in semiconductors with tiny band gaps or semimetals with very small band overlaps. Recently, experimental EI fingerprints have been reported for low-dimensional transition metal dichalcogenide structures exhibiting small band gaps 13,14 ; however, stabilization of a bulk EI state remains elusive. Owing to its semiconductor nature, absence of structural transformations and marked quantum character, ultracompressed 4 He appears to be an excellent candidate in which a bulk EI state could emerge and genuine quantum many-body phenomena like high-temperature excitonic superconductivity and BEC-BCS crossover might exist 15,16 . Is possibly solid helium a bulk EI in the TPa regime? And secondly, the substantial electron-phonon coupling and semimetal Fermi surface previously disclosed in solid helium suggest the possibility of superconductivity in this quantum crystal upon band-gap closure. Is metallic helium a superconductor? If so, what are the underlying physical mechanisms and corresponding critical temperature? Besides their fundamental interest, answering to these questions may have major consequences in the fields of planetary science and astrophysics since this new knowledge could improve our understanding of the thermal and chemical evolution of low-mass WDs 1,2 .
In this study, we employ theoretical first-principles approaches based on DFT and many-body perturbation GW calculations to advance knowledge on the electronic, elastic and superconductor properties of solid helium in the TPa regime. Our main finding is an unprecedented bulk excitonic insulator to superconductor phase transition driven by pressure in which the superconductor state can reach a critical temperature of ≈ 70 K under a compression of 100 TPa. It is worth noting that an exhaustive random sampling of the structural space of solid helium was performed at P = 100 TPa (AIRSS 17,18 ), with the finding that the hcp phase imperturbably remains the ground state (Methods). (Despite the fact that in the TPa regime the c/a ratio of this phase significantly departs from the ideal value of 1.63, we keep labeling it as hcp since its crystal symmetry does not vary.)

Results and discussion
We started by benchmarking different families of DFT functionals (i.e., semi-local, van der Waals corrected and hybrid) 11 against the metallization pressure of solid helium calculated with DMC methods, which amounts to 25.7 TPa (Fig. 1a). In all the analyzed cases, the band gap decreases almost linearly under increasing pressure due to the steady enhancement of electronic delocalization among neighboring atoms (Fig. 1b). The semi-local PBE functional predicts a metallization pressure of 17 TPa, in consistent agreement with previous computational studies 7,8 . Meanwhile, the hybrid functional B3LYP performs the best in comparison to the DMC benchmark by providing a metallization pressure of ≈ 23 TPa ( Fig. 1a and Fig. S1). Van der Waals corrections turn out to be practically negligible in the TPa regime (e.g., the two PBE and PBE-D3 curves practically overlap each with the other) due to the dominant role of interatomic repulsive interactions at short distances 19 . Based on these results, we adopted the hybrid functional B3LYP for our subsequent analysis of the electronic band structure of solid helium.
Unlike atomic hydrogen, hcp 4 He presents an indirect band gap with the VBM located at the reciprocal point Λ and the CBM at the center of the Brillouin zone (Γ point, Fig. 1c), in consistent agreement  4 He. a Electronic band gap, E g , expressed as a function of pressure and calculated with different DFT functionals. Negative E g values indicate overlapping between the VBM and CBM levels. The grey region indicates the stability range of metallic helium as calculated with QMC methods 7 . b Electronic localization function (ELF isosurface = 0.8, yellow) of solid helium at 12 and 36 TPa in a red-green-blue color scale with red denoting high electronic density and blue low electronic density. c Lowest conduction (red) and highest valence (blue) bands expressed as a function of reciprocal wave vector in the k z = 0 plane; the grey and transparent plane represents the Fermi surface.
with previous DFT calculations and experiments 8,10 . It is noted that the direct band gap at the Λ point actually increases under compression (Fig. S2). Interestingly, when the energy gap between the VBM and CBM levels disappears a semimetal state characterized by an almost negligible density of states at the Fermi level emerges due to the fact that no additional electronic bands cross the Fermi surface ( Fig. 2a-c). At 1 TPa, the VBM consists exclusively of s-like orbitals while the CBM exhibits full p-like character (Fig. 2a). Upon further compression, the VBM presents increasingly larger hybridization between s and p-like orbitals while the CBM conserves its pure p-like character (Fig. 2b). At pressures higher than 23 TPa, electrons from the VBM at the Λ point are transferred to the CBM at Γ in order to lower their energy, thus rendering a p-type semimetal system (Fig. 2c).
The continuous pressure-driven closure of the band gap and subsequent stabilization of a semimetal state in hcp 4 He, suggest the possibility of spontaneous formation of excitons with finite momentum |q| = Λ-Γ at low temperatures. An exciton is a bound state formed by an excited electron (e − ) in the conduction band and a hole (h + ) in the valence band that interact through attractive Coulomb forces. In narrow-gap semiconductors, a sufficient condition for the spontaneous formation of excitons is that the corresponding binding energy, E bind , is larger in absolute value than the band gap since then the total energy of the system can be lowered by promoting electrons to the conduction band in the absence of optical excitations 16,20 . We computed the binding energy of an exciton in ultra-compressed hcp 4 He by relying on the Wannier-Mott model ("Methods") since the dielectric constant of solid helium in the TPa regime is relatively high (ϵ r > 5, Fig. S3) and consequently electric field screening effects are large 21 .
Our excitonic binding energy results obtained with the hybrid B3LYP functional and expressed as a function of pressure are shown in Fig. 2d. It was found that the sufficient condition for spontaneous formation of excitons, namely, |E bind | > E g , was fulfilled over a wide pressure interval of approximately 0.2 TPa prior to metallization (cyan region, in which |E bind | B3LY P = 0.19 eV). In view of this result, we performed many-body perturbation theory calculations within the GW approximation to explicitly and more accurately determine quasiparticle excitations in ultracompressed 4 He ("Methods") 22 . As it is shown in the inset of Fig. 2d, GW calculations provided a much larger excitonic binding energy than calculated with the Wannier-Mott model and hybrid DFT functionals, namely, |E bind | GW = 0.34 eV. Moreover, the estimated pressure interval in which excitons can spontaneously form noticeably increased up to 0.45 TPa (orange region). Therefore, based in our hybrid DFT and many-body perturbation GW calculations we may conclude that on the verge of metallization hcp 4 He is a bulk excitonic insulator (EI). The same conclusion was reached when considering alternative structural phases for ultra-compressed solid helium along with semi-local DFT functionals (Fig. S4).
The emergence of a bulk EI state is expected to be accompanied by strong lattice distortions and instabilities due to arising electronphonon interactions 23,24 . We computed the phonon spectrum of hcp 4 He at different pressures using the semi-local PBE and PBEsol functionals (since phonon calculations at this level of theory are feasible), as shown in Figs. 3a and Fig. S5. (The hybrid B3LYP functional certainly provides very similar phonon spectrum results than obtained with semi-local DFT functionals when employing small supercells, Fig. S6.) Reassuringly, a distinct phonon softening appears at the reciprocal lattice point Λ between 15 and 20 TPa, that is, when semi-local DFT functionals predict that solid helium becomes a metal (Fig. S7). Interestingly, above 30 TPa additional phonon softenings emerge along the K-Γ and M-K reciprocal space directions; we found that around this pressure the energy gap between the CBM (located at Γ) and secondary VBM (located at Ω, Figs. 2c and Fig. S7) vanished. Thus, in addition to validating our prediction for the stabilization of a bulk EI state, these findings corroborate the strong coupling between electrons and lattice vibrations previously disclosed in ultra-compressed solid helium 8 . It is worth noting that quantum anharmonic effects were assessed for 4 He with the stochastic self-consistent harmonic approximation (SSCHA) method [25][26][27][28] and found to be of little relevance in the TPa regime ("Methods" and Fig. S8). Besides lattice dynamics, the elastic, structural and thermodynamic properties of solid helium were also found to be influenced by the pressure-driven EI to metal phase transition (Fig. S9). In the insulating phase, the elastic constants of hcp 4 He display a practically linear dependence on pressure whereas in the metallic phase they depart from this behavior and in some cases do not even display a monotonic increase under compression (e.g., C 13 ). Similar effects were also observed for the bulk and shear moduli, sound velocities, Debye temperature and heat capacity (Fig. S9). Regarding the structural features, it was found that the pressure evolution of the hcp lattice parameter ratio c/a drastically changes when the metallic phase is stabilized (Fig. S10). It is worth noting, however, that equivalent alterations were not perceptible in the zero-temperature equation of state (Fig. S11). Thus, taking into account these unanticipated structural, elastic and thermodynamic effects could have important consequences on current modeling of astrophysical bodies, in particular, of low-mass and helium-rich WDs 1,2 .
Motivated by the findings described above, we explored the superconducting properties of ultra-compressed hcp 4 He (using the computationally feasible PBE functional). Accurate electron-phonon coupling (EPC) calculations were carried out with the techniques outlined in the Methods section, which essentially involve the Bardeen-Cooper-Schrieffer (BCS) theory of superconductivity and modified Allen-Dynes formula 29 . Figure 3b shows the Eliashberg spectral function, α 2 F, estimated at different pressures, from which the corresponding average EPC strength, λ, can be straightforwardly computed (Methods and Supplementary Information). At 20 TPa, this function exhibits two appreciable peaks: the most prominent appearing at low frequencies stems from the lowest-energy phonon band at the wave-vector q = (0.25, 0.43, 0), which is very close to the reciprocal space point Λ associated with phonon softening; the other peak emerging at higher frequencies stems from the second and third lowest-energy phonon bands at Γ. The EPC strength associated to these phonon modes respectively amount to ≈140 and 40 (Tables S1-S2), which are extremely high as compared to the average λ of the high-temperature superconductors H 3 S 30 and LaH 10 31 (i.e., 1-10). Such giant λ qν values result from huge phonon linewidths and the minute density of electronic states at the Fermi level ("Methods", Supplementary Information and Fig. S7). However, since the number of phonon modes that appreciably contribute to α 2 F (or, equivalently, to λ) is quite reduced, the superconducting temperature estimated at 20 TPa is relatively low, namely, T c = 17 K (Fig. 3c and Supplementary  Information).
Interestingly, upon further compression, when the energy overlap between the conduction and valence bands is enhanced (Fig. 2c), additional peaks appear in the Eliashberg spectral function that noticeably contribute to the average EPC strength, thus raising the superconducting critical temperature. For instance, at a pressure of 50 TPa, when multiple phonon softenings and α 2 F local maxima are observed (Fig. 3a, b), we estimated a substantial superconducting critical temperature of 47 K (Fig. 3c and Supplementary Information). Under higher compression the superconducting critical temperature steadily increases, reaching a peak value of ≈70 K at the maximum pressure of 100 TPa considered in our calculations (Fig. 3c). It is worth noting that within the pressure interval 20 ≤ P ≤ 30 TPa both T c and λ noticeably decrease; this transient effect is due to a dominant Pinduced surge in the Fermi density of electronic states that drastically reduces the α 2 F peaks (Fig. 3b, "Methods" and Table S3).
An analogous EPC strength parameter and T c analysis was carried out for hcp xenon (Fig. S12) since this material is isoelectronic to solid helium and becomes metallic at experimentally accessible pressures of the order of 0.1 TPa. An EI to metal phase transition similar to that disclosed in ultra-compressed 4 He was also found for hcp Xe at 0.14 TPa. A noticeable phonon softening appeared at a higher pressure of 0.19 TPa, coinciding with the closure of a secondary band gap involving a s-like dominant CBM and p-like dominant VBM (i.e., of the same character than the primary band gap in solid helium). The EPC strength and superconducting critical temperature estimated for hcp Xe at 0.14 TPa are 0.75 and ≈10 K, respectively. Thus, bulk Xe seems to be a good candidate material in which to experimentally search for analogs of some of the key theoretical findings revealed in this work for ultra-compressed hcp 4 He.
As it was mentioned previously, EPC effects on the atomic positions significantly increase the zero-temperature metallization pressure of solid helium by several TPa 8 . Consequently, a well-founded question is: how such type of EPC effects may affect the excitonic insulator and T c results presented in this study? To quantitatively answer this question (although approximately, due to the lack of computationally affordable approaches), we proceeded as described next. First, to evaluate the influence of EPC effects on the excitonic binding energy, we generated a number of supercell configurations for a pressure of 18 TPa in which the atoms were displaced along representative thermal lines determined with the method described in work 32 . Subsequently, we calculated E bind for those configurations by using the Wannier-Mott model, along with E g , and computed their average value ("Methods" and Supplementary Information). It was found that EPC effects (1) significantly open the band gap, in agreement with previous DFT calculations 8 , and (2) enhance the excitonic binding energy very substantially by >40% ("Methods" and Supplementary Information). And second, we carried out a perturbative T c calculation in which we combined the phonon spectrum obtained at ≈30 TPa, that is, the pressure at which solid helium becomes metallic when taking into consideration EPC effects, and the electronic band structure obtained at ≈20 TPa by neglecting EPC effects, which renders an incipient metallic state ("Methods" and Supplementary Information). By proceeding so, we obtained a corrected T c for a nominal pressure of ≈30 TPa that was six times higher than the corresponding value obtained when disregarding dynamical effects. Therefore, we may conclude that EPC on the atomic positions, in addition to displace the metallization of solid 4 He towards higher compression, is also likely to enhance the excitonic binding energies and superconductor transition temperatures reported in this study. Figure 3d shows a sketch of the possible phase diagram of solid helium at pressures and temperatures that are relevant to astrophysical studies. At sufficiently high pressures and low temperatures, a bulk EI state is stabilized. Whether in such a state the spontaneously created electron-hole bound pairs form excitonic Bose-Einstein condensates or exhibit excitonic superconductivity close to zero temperature 15,16 , is a matter that we cannot resolve with the DFT-based methods employed in this work (hence the questions marks in the figure). Upon further compression, hcp 4 He becomes a superconductor with a critical temperature that increases under pressure (made the exception of a small pressure interval following metallization, which has been neglected in the figure). At high enough temperatures, superconductor solid helium transforms into a p-type semimetal. These electronic phase transitions significantly impact the structural, elastic, thermodynamic and transport properties of hcp 4 He hence should be taken into consideration in advanced evolutionary models of stellar bodies like white dwarfs.
In conclusion, we have presented a comprehensive first-principles computational study of the physical properties of solid helium in the TPa regime, putting special emphasis on its electronic band-structure features. It was found that over a broad pressure range preceding metallization hcp 4 He becomes a bulk excitonic insulator in which electrostatically bound electron-hole pairs can form spontaneously. This bulk excitonic insulator state could host genuine quantum many-body phenomena like high-temperature excitonic superconductivity and excitonic BEC-BCS crossover, although additional advanced studies are necessary to fully assess these hypotheses. Upon band-gap closure, solid helium transitions into a superconductor state that possesses a critical temperature of the order of 10 1 -10 2 K, depending on compression. This pressure-induced EI to superconductor phase transition is accompanied by several elastic and structural anomalies. Thus, our theoretical findings besides conveying great fundamental interest are also of great relevance to the physics of celestial bodies, in particular, of low-mass WDs mostly containing metallic helium. Furthermore, it is argued that some analogs of the key theoretical findings revealed here for ultra-compressed helium could be experimentally observed in solid xenon.

First-principles calculations outline
Density functional theory (DFT) calculations were performed with the Vienna ab initio simulation package (VASP) 33 . The projector augmentedwave (PAW) method 34 was employed and the 1s 2 electrons in the He atoms were treated as valence. Different families of DFT functionals were tested among which we highlight the semi-local Perdew-Burke-Ernzerhof (PBE) 35 and revised PBE for solids (PBEsol) 36 , van der Waals corrected DFT-D3 37 , non-local dispersion corrected vdW-optB88 38 , vdW-DF-cx 39 , rev-vdW-DF2 40 , and the hybrid HSE06 41 , B3LYP 42 and HSEsol 43 . Different functionals provide analogous equations of state (Fig. S11). A plane wave energy cutoff of 1500 eV was employed along with dense Monkhorst-Pack k-point sampling grids of resolution 2π × 0.025 Å (Fig. S13). The energy and atomic forces in the structural relaxations were converged to within 10 −6 eV and 0.002 eV/Å, respectively. For validation purposes, we compared our band gap results obtained with the PBE functional as implemented in the VASP code with a full-potential (linearized) augmented plane-wave method as implemented in the WIEN2k code 44 (Fig. S14). Phonon calculations were performed with the small displacement method and the PHONOPY code 45 by employing large supercells of 4 × 4 × 4. A phonon spectrum test was performed for the hybrid B3LYP functional considering a relatively small 2 × 2 × 2 supercell (due to computational limitations). We found that the phonon spectra obtained with the hybrid B3LYP and semi-local PBE functionals presented only minor differences around the Γ point (Fig. S6).
The excitonic binding energy was estimated with the Wannier-Mott formula: where m u = (m e ⋅ m h )/(m e + m h ). In the equation above, m e and m h are the effective mass of the electron at the bottom of the conduction band and the hole at the top of the valence band, respectively. Ry represents the Rydberg constant (=13.6 eV), m 0 the rest mass of the electron, and ϵ r the dielectric constant of the system as referred to vacuum. The electron and hole effective masses were computed like the inverse of the second derivative of the conduction and valence band energies with respect to crystal momentum module, |k|, along the reciprocal space path Λ-Γ. The Wannier-Mott formula is a good approximation for the excitonic binding energy of materials possessing high dielectric constants 21 , which is the case of hcp 4 He in the TPa regime (Fig. S3). The elastic tensor was determined at zero temperature by performing six finite lattice distortions and four atomic displacements of 0.01 Å along each Cartesian direction. The adiabatic bulk modulus, K, and shear modulus, G, were obtained by computing the Voigt-Reuss-Hill averages from the elastic tensor. The longitudinal and transverse sound velocities were calculated with the formulas v p = ðK + 4 3 GÞ=ρ and v s = G=ρ Â Ã 1=2 , respectively, where ρ represents the atomic density of the system.

Crystal structure prediction analysis
The ab initio random structure searching (AIRSS) package 17,18 was used to perform crystal structure searches for solid 4 He. The first-principles DFT code CASTEP 46 was employed to perform the underlying electronic structure calculations based on the PBE functional 35 . The searches were performed at 100 TPa, producing approximately 1000 relaxed structures and considering a total of 12 atoms in the simulation cell. The energy cutoff was set to 1000 eV and a specially designed hard on-the-fly potential was employed for the calculations. Under these conditions, it was found that the hexagonal P6 3 /mmc (hcp) phase remained the ground state followed by a rhombohedral R3m phase with a higher relative energy of 0.131 eV/atom. Other energetically competitive structures were an hexagonal P6m 2 (0.140 eV/atom) and a cubic Im3m (0.170 eV/atom) phase.

SSCHA calculations
Quantum anharmonic effects were assessed with the stochastic selfconsistent harmonic approximation (SSCHA) method [25][26][27][28] . All the SSCHA calculations were evaluated at the pressure of 35 TPa and 0 K, conditions at which excellent convergence of the SSCHA minimization has been verified by an extra population including 800 supercell configurations. SSCHA calculations were performed with a 6 × 6 × 3 supercell including 216 atoms, which yields the dynamical matrices on a commensurate q-mesh of 6 × 6 × 3. The trial harmonic dynamical matrices used for initializing the free energy were obtained from the DFPT method as implemented in the Quantum Espresso (QE) code in the corresponding commensurate q-mesh. In the selfconsistent calculations of the supercells, we used the same cutoff energy as the electron-phonon coupling calculations for the primitive cell, but the k-mesh was reduced accordingly and was tested for convergence. In the SSCHA iterations, except the first four populations in which only internal coordinates were optimized to speed up the minimization, the free energy in other populations was minimized with respect to all degrees of freedom of the crystal structure including the internal coordinates and the lattice cell parameters.

Many-body perturbation theory calculations within the GW approximation
The excitonic binding energy was also estimated by means of highly accurate many-body perturbation theory calculations within the GW approximation 47 performed with the Yambo code 48 . For this, we employed the generalized gradient approximation (GGA) as parameterized by PBE together with a planewave basis set and norm-conserving pseudopotential. The kinetic energy cutoff for the wave functions was set to 600 Ry. The Brillouin zone was sampled with a 64 × 64 × 32 k-mesh. Many-body quasi-particle energies were obtained within the GW approximation 49 by considering 8 G 0 W 0 iterations, and the dynamic dielectric function was obtained with the full-frequency method with up to 800 integration points. The exciton energies were calculated by solving the Bethe-Salpeter equation 22 with the Tamm-Dancoff approximation 50 . The static screening in the direct term was calculated within the random-phase approximation with the inclusion of local field effects. We used 2 valence and 3 conduction bands to solve the Bethe-Salpeter equation matrix. For the GW band-structure calculations, we sampled the Brillouin zone with a 16 × 16 × 8k-point grid. A kinetic energy cutoff of 90 Ry was used for the evaluation of the exchange part of the selfenergy and of 150 Ry for the dielectric screening matrix size. About one hundred unoccupied bands were used to build the polarizability and integrate the self-energy. The convergence test on the number of frequency points and perturbation iterations are shown in Fig. S15. The exciton energies were mapped along high symmetry paths for different pressures, as it is shown in Fig. S16 and summarized in Table S4.
Electron-phonon coupling calculations and critical superconducting temperature Electron phonon coupling (EPC) calculations were performed with the Quantum Espresso (QE) code 51,52 by using ultrasoft pseudopotentials (i.e., semi-local PBE 35 ), an energy cutoff of 200 Ry for the kinetic energy and an energy cutoff of 2000 Ry for the charge density (convergence tests are shown in Table S1). The equation of state of hcp 4 He computed with the VASP and QE codes show very good agreement, as it is illustrated in Fig. S10. The electron-phonon matrix elements were calculated in a 16 × 16 × 8 q-point grid with density functional perturbation theory (DFPT) 53 . We adopted a dense and shifted k-point mesh of 80 × 80 × 40 to increase the convergence in the self-consistent calculations. For the EPC calculations, we further increased the kpoint mesh up to 192 × 192 × 96 (convergence tests are shown in Table S2) and to ensure k-point sampling convergence we employed the Methfessel-Paxton scheme with a smearing width of 0.02 Ry. The Dirac deltas on the band energies were substituted by Gaussian functions with a broadening of 0.002 Ry, which are necessary for the calculation of the EPC strength parameter λ (convergence tests are shown in Fig. S17). Convergence tests on the q-point grid sampling are also presented in Table S5. For further validation purposes, we also performed calculations with a hardcore pseudopotential involving a realspace cutoff of r c = 0.37a 0 54 and compared the results with those obtained with ultrasoft pseudopotentials, as it is shown in Fig. S18.
The Eliashberg spectral function, α 2 F(ω), accounts for the coupling between phonons and electrons in the Fermi surface like: where N(E F ) is the density of states at the Fermi level (per unit cell), γ qν the linewidth of the phonon mode ν at the wave vector q, and N qν the total number of qν points in the sum. The critical superconductor temperature, T c , was estimated with three different formulas: the McMillan formula 55 T mAD c = ð1:0061 where μ * is the Coulomb pseudopotential, for which we selected values within the widely accepted range of 0.10-0.13, and the parameters f 1 and f 2 are defined like: with ω 2 = hω 2 i 1=2 , Λ 1 = 2:46 1 + 3:8μ * À Á and Λ 2 = 1:82ð1 + 6:3μ * Þð ω 2 =ω log Þ. Meanwhile, the logarithmic average phonon frequency, ω log , is defined like: and the EPC strength, λ, is proportional to the first inverse momentum of the spectral function, namely: where λ qν = γ qν π_NðE F Þω 2 qν : ð9Þ The λ qν parameter in the equation above corresponds to the EPC strength of the phonon mode at wave vector q and phonon branch ν. All electron-phonon coupling and T c results are listed in Table S3.
To evaluate the changes to the excitonic binding energy arising from electron-phonon coupling, we first calculated the phonon eigenvalues and eigenvectors on a 4 × 4 × 4 q-point grid using the finite displacement method in conjunction with nondiagonal supercells 57 . We then used those phonons to generate atomic configurations in a 4 × 4 × 4 supercell in which the atoms were displaced along thermal lines 32 , providing representative configurations that the atoms adopt at zero temperature. We then calculated the average of the electronic properties of the system over these configurations to obtain an estimate of the changes to the band gap and excitonic binding energy driven by quantum fluctuations of the helium nuclei, as shown in Table S6. Due to the high phonon energies in solid helium, thermal fluctuations only become important at temperatures higher than those of interest in this work 8 .
We also qualitatively evaluated the effects of electron-phonon coupling on the predicted superconductor transition temperature. For this test, we combined the electronic band structure obtained at 20 TPa (i.e., when the band gap closes in the absence EPC effects) with the phonon spectrum obtained at 30 TPa (i.e., when the band gap closes when considering EPC effects) to calculate the resulting Eliashberg function (Fig. S19). (These α 2 F calculations were performed at the PBE level, just like in the previous cases.) It was found that when EPC effects were approximately considered in such a manner, the predicted T c , increased from 1 K (without corrections, Table S3) to 6 K.

Data availability
The data that support the findings of this study are available upon request from C.L. and the corresponding authors.