Introduction

A Weyl (semi)metal is a new topological state of matter that hosts the condensed matter equivalent of relativistic Weyl fermions. Weyl fermions exist as low-energy electronic excitations at Weyl nodes in three-dimensional momentum space, producing exotic physical properties such as unique surface Fermi arcs and negative magnetoresistance.1,2,3,4 The topological Weyl state can be realized by breaking either time-reversal or lattice inversion symmetry. A candidate topological Weyl semimetal is the quasi two-dimensional transition metal dichalcogenide MoTe2.5,6,7,8 The transition to the non-trivial topologically protected crystal state occurs upon cooling from the high temperature 1T′ monoclinic phase to the low temperature orthorhombic Td phase. The transition is driven by c-axis layer stacking order around 250 K.9,10,11,12,13 Upon cooling to the non-centrosymmetric Td phase, Weyl quasiparticles are expected at characteristic electron and hole band crossings in momentum space.

MoTe2 is a candidate TSC in the orthorhombic phase at ambient pressure.14 In the superconducting state, Fermi arcs are proposed to still exist even though the Weyl nodes are completely gapped out by the superconducting gap, and the onset of superconductivity can generate a quantum anomaly analogous to the parity anomaly in quantum electrodynamics.15,16 The application of pressure enhances the superconducting transition temperature from 0.25 to about 8 K14 and extends the superconducting phase over a wide pressure range. Recent studies showed that while in the superconducting state, a phase transition occurs from the orthorhombic Td to the monoclinic 1T′ phase under pressure17 with a large hysteresis region and coexistence of the two phases. A high-pressure muon-spin rotation study of superconductivity in MoTe2 suggested that a topologically non-trivial s+− state likely exists in MoTe2 at pressures up to 1.9 GPa.18 Key components to a TSC state are the crystal structure and band topology under pressure. Existing reports stop short of exploring the symmetry in the superconducting phase under pressure and no calculations exist on the electronic band structure. Using single crystal neutron diffraction, the pressure–temperature phase diagram is mapped out, and combined with band structure calculations, we elucidate the effects of pressure on the electronic band structure topology. The results provide evidence for a topological state that persists under pressure, consistent with the muon work of ref. 18.

MoTe2 can crystallize into several different structures including the hexagonal 2H, the monoclinic 1T′ and the low temperature orthorhombic Td phases.9,10,11 The hexagonal 2H phase (with space group P63/mmc), a semiconductor with an indirect band gap of 0.88 eV, is stable below 900 °C.19,20 The monoclinic 1T′ phase is stable above 900 °C but can be stabilized to room temperature by quenching.9,10,20 In the 1T′ phase, Mo is surrounded by Te in an octahedral environment with Mo shifted off-center. Mo atoms form zigzag chains running along the b-axis that distort the Te sheets resulting in a tilted c-axis with angle β ~ 93° (Fig. 1a).9 Upon cooling from the 1T′ and not from the 2H phase, a first order structural transition occurs into the orthorhombic Td phase (Fig. 1b). This transition comes about as a result of layer stacking order along the c-axis in the Td phase.13

Fig. 1
figure 1

Physical properties. a Crystal structure of MoTe2 in the 1T′ phase projected on the ab-plane, with the zigzag chains marked running along the a-axis. b Unit cell of the Td structure projected on the ac plane. c Temperature-dependent superconducting (SC) shielding [zero-field-cooled (ZFC)] and Meissner (FC) fraction data for MoTe2 in H = 20 Oe at 1.2 GPa. d Magnetic hysteresis loop for MoTe2 at 2 K at 1.2 GPa. Arrows indicate sweeping direction. e, f Contour map of the neutron scattering intensity obtained by scanning across the (2 0 l) Bragg peak along the [0 0 1] direction at 0.15 GPa on cooling (e) and warming (f). The peak labels denoted with D1 (main domain) and D2 (second domain) are used to distinguish the peaks coming from two different domains, and subscripts 1T′ and Td are used to distinguish the peaks in the different phases. g Temperature dependence of the q integrated neutron scattering intensity for \((201)_{1{\mathrm{T}}\prime }^{{\mathrm{D1}}}\) peak (red squares), \((201)_{{\mathrm{T}}_{\mathrm{d}}}\) peak (blue circles), and \((201)_{{\mathrm{T}}_{\mathrm{d}} \ast }\) (green circles). Arrows indicate sweeping directions. h, i Contour map of the neutron scattering intensity obtained by scanning across the (2 0 l ) Bragg peak along the [0 0 1] direction at 0.95 GPa on cooling (h) and at 1.2 GPa on cooling (i). j The phase diagram under pressure. The region between the two phases shows coexistence of both phases plus the Td* phase. Inversion symmetry is broken only in the Td phase while both Td* and 1T′ are centrosymmetric. k A two-dimensional scan in the h0l plane outside the pressure cell with half-integer peaks from the Td* phase appearing along L. The dashed white lines mark the integer peaks along L that correspond to both the Td and Td* phases while the half-integer peaks belong to the Td* phase only

Results

Shown in Fig. 1c, d are the plots of the bulk magnetic susceptibility, χ(T), and magnetization, M, under pressure as a function of temperature and magnetic field, respectively. Our crystal becomes diamagnetic below 3.0 K at 1.2 GPa, marking the onset of superconductivity. Bulk superconductivity is evident in the full shielding diamagnetic signal at 2.0 K. Note that χ(T) was measured along the ab-plane in which the demagnetizing factor is negligible. The crystal exhibits a large shielding fraction, ~100% (zero-field-cooled (ZFC) magnetization), at 2 K and 1.2 GPa. Magnetic hysteresis loops were also measured at 2.0 K and 1.2 GPa as shown in Fig. 1d indicating that MoTe2 exhibits a typical type-II superconducting behavior.

The 1T′ to Td phase transition under pressure was explored by single crystal neutron diffraction performed in the 20L scattering plane, focusing on structural signatures most affected by variations in the angle β between the a- and c-axes. The crystal was aligned using the lattice constants and β of the 1T′ symmetry (a = 6.3263 Å, b = 3.469 Å, c = 13.7789 Å, α = 90°, β = 93.75°, γ = 90°). A typical contour map of the scattering intensity as a function of temperature at constant pressure is shown in Fig. 1e at 0.15 GPa. At 300 K in the 1T′ phase, a very strong Bragg peak at (hkl) = (2 0 1) is observed along with two weak secondary peaks, at (0.4 0 1) and (1.4 0 1). The Bragg peaks are from the two possible domains in the monoclinic symmetry and are marked with superscripts D1 and D2 indicating domain one and two, respectively. One domain (D1) is overwhelmingly more intense than the other and the crystal was aligned using the primary domain. At 0.15 GPa, the 1T′ to Td phase transition is clearly observed upon cooling, with a discontinuous volume change and the transition temperature shifting to a lower value compared to the ambient pressure value (~240 K). At the transition from 1T′ to Td at 220 K, the 1T′ \((201)_{1{\mathrm{T}}\prime }^{{\mathrm{D}}2}\), \((201)_{1{\mathrm{T}}\prime }^{{\mathrm{D1}}}\) and \((202)_{1{\mathrm{T}}^{}}^{{\mathrm{D}}2}\) peaks disappear and the Td peaks, \((201)_{{\mathrm{T}}_{\mathrm{d}}}\) and \((202)_{{\mathrm{T}}_{\mathrm{d}}}\), appear as shown in Fig. 1e. The \((201)_{{\mathrm{T}}_{\mathrm{d}}}\) and \((202)_{{\mathrm{T}}_{\mathrm{d}}}\) peaks appear at differentlpositions corresponding to the change of β from 93.75° (1T′) to 90° (Td). There is only a single domain and one set of peaks in the Td phase. Alignment of the single crystal remained the same as at 300 K during cooling. On warming, a similar set of Bragg peaks is observed and shown in the contour map of Fig. 1f as on cooling but with an additional Bragg peak that arises from the new Td* phase. While only one peak at \((203)_{{\mathrm{T}}_{\mathrm{d}} \ast }\) is present in this narrow range of the 20 L scan inside the pressure cell, other Td* peaks have been observed along L. This is shown in the 2D plot of the H0L scattering plane from a single crystal measurement at 300 K outside the pressure cell (Fig. 1k). Td* peaks appear at half-integer L-points and are prominent along 20 L and 30 L. No Td* peaks are visible on cooling which is consistent with earlier single crystal neutron diffraction measurements (see ref. 21). A thermal hysteresis is evident with the transition temperature shifting to 250 K on warming, which is higher than that on cooling, consistent with the thermal hysteresis observed at ambient pressure and in transport measurements.21 The temperature dependence of the integrated intensities for \((201)_{1{\mathrm{T}}\prime }^{{\mathrm{D1}}}\), \((201)_{{\mathrm{T}}_{\mathrm{d}}}\) and \((201)_{{\mathrm{T}}_{\mathrm{d}} \ast }\) Bragg peaks are shown in Fig. 1g, where coexistence of the two phases appears in the hysteresis region while the Td* phase only exists in the coexistence region. The structure parameters for the Td* phase obtained from the data outside the pressure cell are provided in Table 1. The Td* phase grows at the expense of the Td phase. At no point in pressure or temperature do we observed the Td* phase becoming the only phase, indicating that Td and Td* coexist. Td* disappears upon further warming into the 1T′ phase. At 0.15 GPa and 230 K data, all three phases coexist (Fig. 1g).

Table 1 Crystal unit cell parameters

As the pressure increases to higher values, a precipitous decrease of the 1T′ to Td transition temperature occurs, eventually leading to the disappearance of the Td Bragg peaks by about 1.2 GPa. The intensity plot of Fig. 1h at 0.95 GPa clearly shows the phase transition temperature going down until its complete disappearance by 1.2 GPa as shown in Fig. 1i. At 1.2 GPa, only 1T′ Bragg peaks are observed in the temperature range from 5 to 300 K, regardless of the thermal cycling. Note that the two-domain feature of the 1T′ phase is observed at all pressure points. The single crystal results were additionally confirmed by neutron powder diffraction measurements under pressure that reached down to 1.5 K. The phase diagram under pressure is shown in Fig. 1j. Shown on the diagram are three distinct regions: the high temperature monoclinic region, the low temperature orthorhombic region and the coexistence region between the two that is also host to the new Td* phase. On warming, the Td* phase appears up to about 1 GPa; unit cell doubling occurs in the coexistence region of the two phases that disappears upon applying pressure greater than 1 GPa. This is shaded in green on the phase diagram. Above 1 GPa, a broad transition from the Td to 1T′ with apparent phase coexistence is observed. Coupled with the transition to the monoclinic phase under pressure is suppression of the magnetoresistance22 which may be linked to the disappearance of the Weyl nodes.

Shown in Fig. 2 is a comparison of the Fermi surface projected on the kxky momentum plane in the orthorhombic Td phase, calculated at two pressure points, 0 GPa (upper panels) and 0.15 GPa (lower panels). At 0 GPa, the two Weyl points indicated by dots in the kxky plane at kz = 0 and at E − EF = 23 and 52.5 meV appear at the intersection of electron and hole pockets as indicated in both figures and also in the band structure diagram on the right where the crossing points in the dispersion are marked by the letters A and B. These Weyl nodes are the same as those reported in ref. 21. Increasing the pressure to 0.15 GPa changes the c-axis lattice constant by about 0.5% (with negligible change in the a-axis lattice constant). At 0.15 GPa, new nodes appear at E − EF = 18 and 58 meV above the Fermi level, on the same band as the nodes at E − EF = 23 and 52.5 meV at 0 GPa, but at kz ≠ 0. At 0.15 GPa, the superconducting transition temperature increases as reported in ref. 14 while magnetoresistance goes down.22 The same nodes could not be traced in the band dispersions from calculations performed at higher pressures. At pressure values exceeding 1.2 GPa, the 1T′ symmetry is centrosymmetric and no Weyl nodes are expected. Ambient pressure calculations confirmed that the 1T′ with the P21/m space group remains topologically trivial.13 Thus, with warming from the Td to Td* and then to the 1T′ phase, the type-II Weyl nodes are destroyed.13

Fig. 2
figure 2

Electronic Structure: The top panels are calculations performed in the Td phase at 0 GPa. The lower panels are calculations performed at 0.15 GPa in the Td phase. The Fermi surface cuts shown on the left panels reveal the Weyl nodes at the intersections of electron and hole pockets. These are marked by A and B. The right panels are band structure plots showing the bands crossing at A and B

Discussion

Analysis of the Berry curvature of the new nodes that appear in the Td phase from the band structure calculations performed at 0.15 GPa shown in Fig. 2 indicates that indeed these new crossings are not trivial. This is demonstrated in Fig. 3 which is a comparison of the energy isosurfaces and dispersion relations of the Weyl nodes at P = 0 and at 0.15 GPa in the Td phase. The Fermi surface morphology changes under pressure as seen from Fig. 2 while the Weyl nodes, observed closed to the Γ point and indicated by A and B in the plots of Fig. 2, change energy. Calculations of the energy isosurface at the nodes based on k·p perturbation theory suggest that the Weyl nodes are still type-II as shown in the line plots of Fig. 3 for the two Weyl nodes in the Td phase at 0 and 0.15 GPa at the corresponding k-points, but become significantly tilted under pressure, especially for the second Weyl node at 18 meV. The values of E, kx, ky, and kz for the two Weyl points are provided in Table 2. These results suggest that with the enhancement of superconductivity, the phase remains topologically non-trivial. Moreover, changes in the Berry curvature may be linked to the decrease in magnetoresistance under pressure.22

Fig. 3
figure 3

k.p calculations: Energy isosurfaces and dispersion relations of the Weyl nodes in Td phase are shown in the figure: ah Td phase at 0.0 GPa, ip Td phase at 0.15 GPa. The isosurfaces have the same energy as the corresponding Weyl nodes, and the dispersion plots are along the three principle axes of the matrix M. The existence of such energy isosurfaces and dispersion relations indicates that the Weyl nodes are type-II

Table 2 Electronic structure parameters

Last, we comment on the prospect of topological superconductivity in the material. While superconductivity is observed in MoTe2 in the Td and 1T′ phases under a range of pressures, topological superconductivity is not guaranteed unless the following two conditions are satisfied: In (i), the Weyl fermions have to either be responsible or take part in the superconducting pairing. And in (ii), assuming the pairing order is s-wave and does not break time-reversal symmetry spontaneously, the sign of the pairing order parameter must be non-uniform. Reference 23 shows that the integral topological index N, which counts the (net chiral) number of surface Majorana cones, is related to the combination of signs \(\mathop {\sum}\nolimits_j {c_j} {\mathrm{sgn}}({\mathrm{\Delta }}_j)/2\), where cj is the chirality of the jth Weyl fermion and sgn(Δj) is the sign of its pairing order parameter in the Bogoliubov-de Gennes Hamiltonian. Trivial Fermi surfaces that have trivial Chern numbers do not affect the topology of the superconducting state assuming the pairing order is BCS s-wave. First, the type-II nature of the Weyl cones in MoTe2 provides a huge hyperbolic Fermi surface and density of state at the Fermi level, both of which favor condition (i). Second, thanks to the small crystal symmetry group of the material, the Weyl points are not related to each other, other than by time-reversal. In contrast with Weyl (semi)metals with large crystal symmetry group, such as the TaAs family, the pairing order parameters of different Weyl fermions are related by lattice rotation and/or mirror symmetries assuming they are not spontaneously broken. In order of these parameters to take non-uniform value, a strong symmetry breaking field or deformation is needed. On the other hand, for the current small symmetry present in MoTe2 material class, there is no such requirement. Consequently, the pairing response to deformations, such as thermal and crystal structural changes, can be different between unrelated Weyl fermions, and may lead to non-uniform pairing signs through phase transitions driven by external tuning parameters. For instance, in Table 2 we show that the tilting behavior of the Weyl points along the principal axes can change as a function of pressure. While conditions (i) and (ii) as well as the protected surface Majorana fermions are yet to be verified, MoTe2 can be a promising material for further investigation due to its Weyl electronic structure’s strong dependence on thermal and structural parameters.

Methods

Meaurements

MoTe2 single crystals were grown by self-flux method. The details of single crystal growth have been described elsewhere.16 The typical size of the obtained single crystals is about 4 × 2 × 1 mm3. Single crystal neutron diffraction under high pressure is combined with electronic band structure calculations to investigate the superconducting transition. Our MoTe2 single crystal exhibits a superconducting transition near 3 K at 1.2 GPa from bulk magnetic susceptibility measurements under pressure. The bulk magnetic susceptibility was measured at high pressure using a vibrating sample magnetometer (VSM) to demonstrate the superconducting transition in our MoTe2 single crystals. The applied pressure was calibrated by measuring pressure cell compression and the pressure dependence of the superconducting transition temperature of a tiny lead piece. High-pressure neutron diffraction experiments were performed on single crystals at HB-1A triple axis spectrometer and using powders at HB-2A powder diffractometer at the High Flux Isotope Reactor (HFIR) of Oak Ridge National Laboratory (ORNL). Self-clamped piston cylinder cells made with CuBe was used as pressure cells for both bulk magnetic susceptibility measurements and neutron diffraction experiments. Daphne 7373 oil was chosen as the pressure transmitting medium. Single crystal neutron diffraction were performed at HB1 triple axis spectrometer with a fixed incident neutron energy of 13.5 meV. The lattice constant measurements of NaCl were used to calibrate the pressure. The single crystal neutron diffraction data shown in Fig. 1k were collected at the CORELLI spectrometer at ORNL.

Calculations

The theoretical DFT calculation was implemented in the Vienna Ab initio Simulation Package (VASP) for the pressurized Td and 1T′ phases of MoTe2. VASP is a plane-wave based projector augmented pseudopotential local density functional method. Its exchange-correlation functional takes the generalized gradient approximation parameterized by Perdew-Burke-Ernzerhof. We turned on the spin–orbit coupling, used the default energy cutoff for the calculation and chose a K-grid of 13*7*3. An NBANDS of 160 is used to allow enough bands for Wannier90 to fit a minimum basis tight binding Hamiltonian involving Mo’s 4d and Te’s 5p orbits with necessary symmetrization. We used Wannier_tools for ab initio data processing, including nodal searching, Chirality calculation, and Fermi surface contour plot. For the Td phases, the lattice constant for 0 GPa is a = 6.2977 Å, b = 3.474 Å, c = 13.8286 Å, the lattice constant for 0.15 GPa is a = 6.3213 Å, b = 3.461 Å, c = 13.7499 Å. For the 1T′ phases, the lattice constant for 0 GPa is a = 6.33 Å, b = 3.469 Å, c = 13.86 Å with β = 93°, the lattice constant for 0.15 GPa is a = 6.3263 Å, b = 3.467 Å, c = 13.7789 Å with β = 93.75°, the lattice constant for 1.2 GPa is a = 6.2954 Å, b = 3.4601 Å, c = 13.4766 Å with β = 94.183°, the lattice constant for 1.8 GPa is a = 6.294 Å, b = 3.456 Å, c = 13.381 Å with β = 94.183°. The coordinates of Mo and Te atoms at 0 GPa are obtained from the refinement of the neutron scattering results while for systems under pressure they are simply rescaled to the new lattice constant obtained under pressure.

In the k.p perturbation theory, the Hamiltonian is expanded up to the first order term near each of the Weyl nodes. In the basis of the two bands that crosses at the Weyl node, the Hamiltonian can be written as a 2 by 2 matrix. A 2 by 2 Hermitian matrix has four independent variables, so for a 3D case the p part of the Hamiltonian will have totally 12 independent variables. With the help of Pauli matrices the k.p Hamiltonian can be expressed in the following form:

$${\mathrm{H}}(\overrightarrow {\mathrm{w}} + {\mathrm{\delta }}\overrightarrow {\mathrm{k}} ) = [{\mathrm{\varepsilon }}(\overrightarrow {\mathrm{w}} ) + {\mathrm{\delta }}\overrightarrow {\mathrm{k}} \cdot \overrightarrow {\mathrm{v}} _0]{\mathrm{\sigma }}^0 + {\mathrm{\delta }}\overrightarrow {\mathrm{k}} \cdot {\mathrm{R}} \cdot \overrightarrow {\mathrm{\sigma }} ^{\mathrm{T}}$$

Here \({\vec{\mathrm w}}\) is the reciprocal vector of a Weyl node, \({\mathrm{\varepsilon }}({\vec{\mathrm w}})\) is the energy at that node, σ0 is the 2 × 2 identity matrix, \({\vec{\mathrm \sigma }}^{\mathrm{T}} = ({\mathrm{\sigma }}^1,{\mathrm{\sigma }}^2,{\mathrm{\sigma }}^3)^{\mathrm{T}}\) is the vector formed by the three Pauli matrices, \({\vec{\mathrm v}}_0\) and R are real 1 by 3 vector and 3 by 3 matrix that contain the 12 independent varibles. The energy dispersion near the Weyl node can then be solved analytically:

$${\mathrm{E}}^ \pm (\overrightarrow {\mathrm{w}} + {\mathrm{\delta }}\overrightarrow {\mathrm{k}} ) = {\mathrm{\varepsilon }}(\overrightarrow {\mathrm{w}} ) + {\mathrm{\delta }}\overrightarrow {\mathrm{k}} \cdot \overrightarrow {\mathrm{v}} _0 \pm \sqrt {{\mathrm{\delta }}\overrightarrow {\mathrm{k}} \cdot {\mathrm{M}} \cdot {\mathrm{\delta }}\overrightarrow {\mathrm{k}} ^{\mathrm{T}}}$$

Where \({\mathrm{M}} = 3{\mathrm{R}} \cdot {\mathrm{R}}^{\mathrm{T}}\) is a real symmetric matrix. The group velocity along any given direction on each of the two crossing bands at each Weyl nodes can be calculated using the DFT method described above. The vector \({\vec{\mathrm v}}_0\) and the matrix M can then be determined from these velocity data.