Superconductivity from the condensation of topological defects in a quantum spin-Hall insulator

The discovery of quantum spin-Hall (QSH) insulators has brought topology to the forefront of condensed matter physics. While a QSH state from spin-orbit coupling can be fully understood in terms of band theory, fascinating many-body effects are expected if it instead results from spontaneous symmetry breaking. Here, we introduce a model of interacting Dirac fermions where a QSH state is dynamically generated. Our tuning parameter further allows us to destabilize the QSH state in favour of a superconducting state through proliferation of charge-2e topological defects. This route to superconductivity put forward by Grover and Senthil is an instance of a deconfined quantum critical point (DQCP). Our model offers the possibility to study DQCPs without a second length scale associated with the reduced symmetry between field theory and lattice realization and, by construction, is amenable to large-scale fermion quantum Monte Carlo simulations.


I. SUPPLEMENTARY NOTES
Observables. Symmetry-broken states are characterised by a local order parameterÔ r,δ , where r denotes a unit cell and δ an orbital within the unit cell. The associated time-displaced correlation functions read For the finite-size scaling analysis, we consider the order parameter the equal-time correlation ratios with |∆q| = 4π √ 3L , the susceptibilities and the corresponding correlation ratios Here, Λ 1 () indicates the largest eigenvalue of the corresponding matrix in orbital space (6 × 6 for spin currents, 2 × 2 for pairing) and the ordering wave vector is at the Γ point. These quantities exhibit the following finite-size scaling behaviour near the critical point: Here, λ c , ν, η, and z are the critical coupling, the correlation length exponent, the anomalous dimension, and the dynamical critical exponent, respectively. The correlation ratios R O and R O χ are both renormalization group (RG) invariant quantities at the critical point and hence provide a simple way to estimate λ c and ν without any knowledge about η. However, the generic corrections-to-scaling exponent ω is not necessarily the same for all four quantities in equation (6). Such corrections generally arise from irrelevant operators of the fixed point and the analytic part of the free energy. If the absolute value of the negative RG dimension is relatively large, the main contribution to ω will come from the background term of the free energy [1]. In this case, This suggests that the susceptibility χ and the corresponding correlation ratio R χ will have smaller scaling corrections than the corresponding equal-time quantities if the effect of the negative RG dimension is small at λ c . We assumed a dynamical critical exponent z = 1 for both the SM-QSH and the QSH-SC transition. This is motivated by the Lorentz invariance of the corresponding field theories [2,3]. Accordingly, in our simulations, we used β = L and thereby fixed L z /β in the above finite-size scaling expressions.
Trotter decomposition used for the QMC simulations. For the finite-temperature QMC simulations underlying this work, imaginary time was discretized with a spacing ∆τ = β/L τ . To ensure hermiticity, the partition function is written as whereĤ t is defined by equation (1) and the interaction (2) was partitioned into local operators H α λ (i) acting on spin component α = x, y, z and on hexagon i. The leading discretization error for the partition function then scales as ∆τ 2 .
The Trotter decomposition in equation (8) breaks the global SU(2) spin rotation symmetry. For example, [Ĥ x λ (i),Ĥ y λ (i)] = 0, so that equation (8) will not be invariant under a global SU(2) rotation. Because SU(2) symmetry breaking is a relevant perturbation for both critical points considered, care has to be taken to ensure that its effects, which scale as ∆τ 2 , remain below the relevant energy scale. An explicit test involves the total spin operator and generator of global SU(2) rotations Here, S r,δ =ĉ † r+δ σc r+δ and δ runs over the positions of atoms in the unit cell at r. Supplementary Fig 1 shows the associated time-displaced spin-spin correlation function. A global SU(2) spin symmetry implies that this quantity is independent of imaginary time. The numerical results are essentially constant in imaginary time if the symmetric Trotter decomposition is used. Therefore, the latter was employed together with ∆τ = 0.2 for all results of this work.

II. SUPPLEMENTARY METHODS
Finite-size scaling analysis by different crossing points. In addition to the main text, where a 'two-size crossing' with sizes L and L + 6 was performed, we provide a cross check based on crossings of L and L + 3.
When considering the L and L + 6 crossings in the main text for the estimation of λ SC c2 (see Fig. 3 of the main text), we are obliged to take into account the L = 6 data. Upon inspection, the fit turns out to be rather bad since χ 2 /DOF = 66.9. On the other hand, if we consider the L and L + 3 crossing points, we can omit the L = 6 data and get a more acceptable χ 2 /DOF = 6.8. As apparent from Supplementary Fig. 3, the extrapolated value of λ SC c2 based on the crossing points of L and L + 3 compares favourably with the analysis in the main text. Supplementary Fig 2a,b show the crossing values of λ and 1/ν at the semimetal-QSH transition, as obtained from the correlation ratio. Supplementary Fig 2 c shows the anomalous dimension η at each crossing point, as well as a fit based on an expression analogous to Eq. (6) of the main text, where r = L+3 L . The results of a similar analysis for the QSH-SC transition are reported in Supplementary Fig. 3.
Consistency of the finite-size scaling analysis. As an independent consistency check on the results from the 'two-size crossing' method used in the main text, we consider in this section a collective fitting of multiple system sizes.

DQCP -QSH
To reduce the number of degrees of freedom in the fit of the anomalous dimension η, a substitution in terms of the scaling form for χ O and R O χ in equation (6) is performed, using the expansion (we ignore the correction-to-scaling term) As shown in table 2, the collective fitting for both order parameters is acceptable for L min ≥ 12. The corresponding values η QSH = 0.194(9) and η SC = 0.279(9) agree well with the values 0.21(5) and 0.22(6) obtained with the 'two-size crossing' approach used in the main text.
We also carried out the collective fitting at the Gross-Neveu critical point. In contrast to the DQCP, this phase transition suffers much less from corrections to scaling (as shown in Fig. 3 of the main text, the crossing points converge quickly). Hence, a fit without the scaling correction term is performed, and the results are shown in table 3. As can be seen, a good χ 2 /DOF is obtained once the L = 6 data set is neglected, and fits for L min = 9, 12, 15 or 18 produce consistent results. Taking λ c1 = 0.01891(2) and 1/ν = 1.17(3) from L min = 9, the results match those of the analysis in the main text. Results for the anomalous dimension at the Gross-Neveu critical point are listed in table 3. The fits yield acceptable χ 2 values for L min ≥ 12. The exponent η = 0.76(1) from L min = 12 also matches the analysis in the main text.
III. SUPPLEMENTARY DATA Single-particle gap and free-energy derivative across the QSH-SC transition. The single-particle gap ∆ sp is obtained from the single-particle Green function where r + δ runs over the two orbitals of the unit cell located at r. The single-particle gap is minimal at the Dirac point K = ( 4π 3 , 0) and is extracted by noting that asymptotically Here, we used β = 36. Supplementary Fig 4a demonstrates that ∆ sp remains nonzero across the QSH-SC transition at λ c2 ≈ 0.033. In order to clarify nature of the QSH-SC transition, we also calculated the first partial derivative of the free energy density with respect to the coupling λ (we use the same notation as in the main text) Supplementary Fig 4b shows ∂F/∂λ for β = L in the vicinity of λ c2 ≈ 0.033. As expected for a continuous transition, we observe no sign of a jump.

IV. SUPPLEMENTARY DISCUSSION
Charged skyrmion defects of the QSH state. In Ref.
[4], it was shown that when a QSH state is generated by spontaneous symmetry breaking, skyrmion defects of the vector order parameter will carry an electric charge of Q e = 2e, leading to the relation where Q is the Pontryagin index that counts the winding of the unit vector order parameter on the sphere.
Here, we substantiate this fact in terms of an explicit calculation for a lattice model. Our starting point is the Hamiltonian where N (x) = (N x (x), N y (x), N z (x)) is a unit vector at position x corresponding to the centre of a hexagon. SinceĤ is invariant under time reversal symmetry,T −1 α ĉ i,↑ c i,↓ T = α ĉ i,↓ −ĉ i,↑ , Kramers' theorem holds and stipulates that all eigenstates are doubly degenerate.
On the honeycomb lattice, the Pontryagin index is defined as with unit vectors a 1 = (1, 0) and a 2 = ( 1 2 , √ 3 2 ). For an arbitrary vector field N (x), Hamiltonian (17) does not preserve particle-hole (P-H) symmetry. For example, defining the P-H transformation asP where A general P-H transformation can be written aŝ For Hamiltonian (17), it yieldsP where and Thus, there is no generic P-H transformation that leaves this Hamiltonian invariant, unless N (x) is varied in an R 2 space (θ and φ can be defined such that sin θ cos φN x + sin θ sin φN y + cos θN z = 0). Since the transformation has a determinant of −1, the generic P-H transformation gives The 'electric charge' refers to the number of occupied states at zero temperature, relative to half filling. The sign change of the Pontryagin index under a P-H transformation provides a natural way of understanding equation (16). Note that in contrast to a skyrmion, a 2D topological defect (such as vortex) has a vanishing Pontryagin index and carries no charge. The argument of charged skyrmions fails when other Dirac mass terms are considered. For example, a system with fluctuations of a three component vector field Yukawa-coupled to the three antiferromagnetic mass terms does not break P-H symmetry. In this case, a skyrmion configuration with nonzero Pontryagin index does not carry electric charge.
We diagonalised Hamiltonian (17) on a honeycomb lattice with L = 36, setting t = 1 and λ = 0.5. Supplementary Fig 5  compares the density of states for a uniform field N (x) and for a 'hedgehog' configuration corresponding to a single skyrmion. On the lattice, the Pontryagin index is not quantised and we obtain Q ≈ −0.989. The system remains gapped when one skyrmion is inserted, see Supplementary Fig. 5b. The breaking of P-H symmetry is also apparent from D(ω) = D(−ω). Simple number counting shows that Compared to the case of uniform polarisation in Supplementary Fig. 5a, an additional charge 2e is generated.
On a system with open boundary conditions, the Pontryagin index is not necessarily quantised and we can investigate how charge is transferred during the insertion of a skyrmion by varying the Pontryagin index from zero to one. Supplementary Fig 6 shows that the total charge is 'pumped' from 0 to 2e and we observe a step function at a non-integer value of Q. This is a consequence of the aforementioned Kramers theorem. As shown in Supplementary Fig. 7, the bulk remains gapped during this process, while the edge stays gapless. Thus, the charge 2e is pumped through the edge under insertion of a skyrmion. . b Time-displaced spin correlation function for the symmetric Trotter decomposition and λ = 0.019, β = L. Reported error bars correspond to standard errors and are smaller than the symbol size.  Figure 3. Deconfined QSH-SC transition. a Estimation of the critical values λ QSH c2 = 0.03317(4) and λ SC c2 = 0.0332(1). b,c Critical exponents 1/ν SC = 2.1(7), 1/ν QSH = 2.0(2), η SC = 0.21(5), and η QSH = 0.19(9) from finite-size scaling of the crossing points of L and L + 3. Reported errors and error bars correspond to standard errors.