Superconductivity near a ferroelectric quantum critical point in ultra low-density Dirac materials

The experimental observation of superconductivity in doped semimetals and semiconductors, where the Fermi energy is comparable or smaller than the phonon frequencies, is not captured by conventional theory. In this paper, we propose a mechanism for superconductivity in ultra low-density Dirac materials based on the proximity to a ferroelectric quantum critical point. We consider both ionic and covalent crystals. We derive a low-energy theory that takes into account the strong Coulomb forces of the polarization modes and the direct coupling between the electrons and the soft phonon modes. We show that in the case of ionic crystals, Coulomb repulsion is strongly screened by the lattice polarization near the critical point. Using a renormalization group analysis in the low-density limit, we demonstrate that the effective electron-electron interaction is dominantly mediated by the transverse phonon mode. We find that the system generically flows towards strong electron-phonon coupling. Hence, we propose a new mechanism to simultaneously produce an attractive interaction and suppress strong Coulomb repulsion, which does not require retardation. We obtain qualitatively similar results for covalent crystals, though the Coulomb screening in this case is much weaker. We then apply our results to study superconductivity in the low-density limit. We find strong enhancement of the transition temperature upon approaching the quantum critical point. Finally, we also discuss scenarios to realize a topological $p$-wave superconducting state in covalent crystals close to the critical point.

A key ingredient in superconductivity is the pairing between electrons. In spite of the strong Coulomb repulsion in free space, at low energy electrons experience an effective attraction and form bound statesthe Cooper pairs. Thus, superconductivity essentially relies on a mechanism that concomitantly reduces the Coulomb repulsion and generates an attractive interaction. In metals, the attractive interaction results from the interchange of longitudinal phonons that couple to the electronic density. This mechanism requires retardation, namely, the crystal vibration must be much slower than the electronic motion. In terms of energy scales, this implies that the Fermi energy is much larger than the Debye frequency. In the intermediate frequency regime, between these two scales, the Coulomb repulsion is logarithmically suppressed, while the phonon interaction is unaffected [1,2]. As a result, the net interaction between electrons may become attractive below the Debye energy.
From this perspective, systems of low carrier concentration, such as doped semimetals and semiconductors, are not expected to exhibit superconductivity. This is because they have a low Fermi energy, which is comparable to, or even smaller than, the typical Debye frequency, and thus does not allow for the dynamical screening of the repulsion. Moreover, the superconducting transition temperature is exponentially sensitive to the density of arXiv:1901.11064v1 [cond-mat.supr-con] 30 Jan 2019 states, which is typically two orders of magnitude smaller in semimetals and semiconductors when comparing with a standard metal.
Surprisingly however, superconductivity in doped semimetals and semiconductors is ubiquitous. It was first discovered in SrTiO 3 [3] and later in many other materials [4]. To the best of our knowledge [5], the lowest density superconductors discovered to date (in reverse order) are Tl doped PbTe [6], Sr doped Bi 2 Se 3 [7], YPtBi [8], SrTiO 3−x [9] and elemental Bismuth [10]. It is noteworthy that, except for SrTiO 3 , all of these materials are either narrow band topological insulators or topological semimetals. The common feature is a nearly touching conduction and valence bands with Dirac-like dispersion.
It is also interesting to notice that SrTIO 3 , PbTe, and SnTe naturally reside close to a paraelectric-ferroelectric phase transition [11,12], which can be tuned in various manners. Recent experiments show that the superconducting transition temperature in SrTiO 3 depends on the distance from the ferroelectric quantum critical point [13][14][15]. Thus, it is natural to study the relation between the ferroelectric quantum critical point and superconductivity.
A variety of theoretical frameworks have been proposed to discuss superconductivity in the limit of low density, including polar phonons [15][16][17][18][19], plasmons [17,[20][21][22], multi-band effects [23,24], soft optical phonons [25], the charge Kondo effect [26], and instantaneous attraction [27,28]. It is particularly important to single out the seminal contribution of the authors of Ref. [16], who pointed out an essential ingredient in any theory of low-density superconductivity: a long-ranged attractive interactions. Only an interaction with a range that is comparable to the interparticle distance can lead to a sufficiently large coupling constant and drive an instability to superconductivity. This is similar to the phenomenon of Wigner crystallization, where long-ranged Coulomb interaction dominates the kinetic energy in the dilute limit rather than the high-density of states limit. Thus, the constraint of long-ranged attraction narrows down the range of viable pairing mechanisms in the lowdensity extreme limit. Such interaction may result from a dynamically screened Coulomb repulsion [16,20], fluctuations of an order parameter close to a quantum critical point [29][30][31][32] or Goldestone mode fluctuations in certain types of spontaneously broken continuous symmetries [33].
It has been theoretically proposed that pairing in SrTiO 3 is mediated by the fluctuations near a ferroelectric (FE) quantum critical point (QCP) [34][35][36]. However, these studies leave a few important questions open. To understand them let us first quickly review some basic facts regarding the FE QCP. This transition is essentially a structural transition where the order parameter is a vector (a lattice distortion) that spontaneously breaks inversion and rotational symmetry in the ordered state. For example, consider the diatomic ionic crystal in Fig. 1. In the ferroelectric phase, the two ions in the unit cell  1. (a) A diatomic cubic ionic crystal in the paraelectric phase. The two ions are represented by blue (+) and red (-) circles (e.g. blue -Pb and red -Te). Unless these two atoms are identical they will have an average charge imbalance. (b) The ferroelectric phase (inversion-breaking) is characterized by a uniform optical phonon displacement vector u = u b −ur.
Due to the charge imbalance between the ions, this also leads to a finite dipolar polarization density P = Qu, where Q is the charge.
are distorted from their cubic Bravais lattice points. Because they have a charge imbalance (ionic), they induce a uniform electric polarization density. Thus, dynamically, the transition is described by a soft optical phonon mode associated with the relative displacement of the two charged ions. Such a phonon mode has three polarizations: one longitudinal-optical (LO) and two transverseoptical (TO) modes. However, in contrast to naive expectation, the dipolar interactions between these distortions prevent the LO mode from becoming soft at the trasition point [37]. Consequently, the soft bosonic modes associated with the dynamiccs of the FE QCP are purely transverse, which are typically only weakly coupled to the Fermi surface. Thus, the main question that remains unanswered by Refs. [34][35][36] is whether these transverse modes can couple strongly to gapless electrons. For more details we refer to a recent comment written by one of us [38].
In this paper, we answer this question. We show, for the particular case of a Dirac material, how the transverse modes couples to gapless electrons in the longwavelength limit (see Sec. II). We then use the renormalziation group (RG) approach to study the tendency of this theory close to the critical point in Section III. We find that the combination of the strong polar dynamics and a nearly touching conduction and valence bands of the Dirac theory leads to strong screening of the Coulomb repulsion alongside with an enhancement of the electronic coupling to the soft transverse phonon mode. Thus, the phonon-mediated attraction by the transverse modes generically overcomes the Coulomb repulsion close to a FE QCP. Interestingly, this flow takes place even at zero electronic density. Therefore, it is distinct from the standard Tolmachev-Anderson-Morel mechanism to obtain attraction [1,2]; it does not require the phonon frequency to be smaller than the Fermi energy. Finally, in Sec. IV, we analyze the possible superconducting instabilities from the interaction mediated by the critical phonon (ferroelectic) mode. We find strong enhancement of the transition temperature T c due to the enhancement of the electron-phonon coupling close to the critical point.
For completeness, in every step of the way, we compare our analysis with the case of covalent crystals (nonionic) where both LO and TO modes are soft at the critical point. In this case, the screening of the Coulomb repulsion is much weaker (only logarithmic). Interestingly, however, the interplay between phonon-mediated attraction and Coulomb repulsion opens a possibility of topological p-wave superconductivity in a certain range of parameters.
Before proceeding to the analysis itself, we would like to note that the case of a Dirac theory is somewhat simpler than the case of a finite Femri surface. This allows us to form a clearer picture of the fate of the transition and puts the calcualtions under control. However, we also draw motivation from a relastic system: the ionic alloy Pb 1−z Sn z Te, which undergoes a FE phase transition at z = z F E ≈ 0.25 [12]. The transition can be described by the spontaneous formation of a relative displacement between the Pb and Te lattices, therefore, it leads to a gapless optical phonon mode (see Fig. 1). When z is further increased above z T ≈ 0.41, the alloy undergoes a second, topological phase transition, between a trivial insulator and a topological crystalline insulator [39]. The topological transition entails gapless Weyl points close to the L-points of the Brillouin zone [40]. When doped with Tl or In atoms, this alloy becomes metallic and superconducting, with the transition temperature T c exhibiting a peak at some intermediate value of z [41].
Thus, it seems doped Pb 1−x Sn x Te is perfectly suitable for our theory as it includes a FE QCP, a small Fermi surface with Dirac dispersion and superconductivity. However, it is important to note that the situation is more complex. In pure PbTe, for example, superconductivity appears only when doped with Tl. Additionally, it has been found that the superconducting state emerges only above a critical density where additional electron pockets become populated [42]. Nonetheless, we find the question of the ferroelectric quantum critical flcutations coupled to gapless fermions an interesting problem, which is definitely relevant to this alloy and possibly relevant to other system.

II. MODEL
We now consider the low energy effective field theory of a Dirac semimetal near the ferroelectric transition. The Euclidean (imaginary time) action is given by the sum where the first three terms describe the dynamics of the fermions ψ, the optical phonon field u and the Coulomb field φ, while the latter three describe their interactions. The Coulomb field φ should be considered as a Hubbard-Stratonovich decomposition of the instantaneous Coulomb interaction. Now we specify these terms in detail.

A. Quadratic terms
The electron term -The electron quadratic term (motived by the model of the PbTe crystal [39]) reads where ψ n is a four-component Dirac spinor, n = 1, . . . , N denotes different fermionic flavors, andψ n ≡ ψ † n γ 0 . Parameters v F , m, and ε F stand for electron velocity, Dirac mass, and Fermi energy, respectively. We use Hermitian gamma matrices {γ 0 , γ 1 , γ 2 , γ 3 } = {σ 1 ⊗ σ 0 , σ 2 ⊗ σ 1 , σ 2 ⊗ σ 2 , σ 2 ⊗ σ 3 } and γ 5 = γ 0 γ 1 γ 2 γ 3 , where σ i are usual Pauli matrices. Notice that here we have assumed an isotropic dispersion by taking the same velocity v F in all directions. The anisotropic case does not modify the main qualitative results of this paper, and therefore we comment on it in Appendix D. For generality, we have assumed a non-zero mass term m and a finite Fermi energy ε F . However, we will neglect them in our RG analysis, assuming that they are much smaller than other relevant energy scales.
We note two important discrete symmetries of Eq. (1): inversion symmetry P and time reversal symmetry T . The action of these symmetries in terms of Dirac matrices is given by P = γ 0 and T = γ 1 γ 3 K, respectively, where K is complex conjugation. The phonon term -Next, we consider the dynamics of the phonon modes, which become soft at the FE phase transition. To have an intuitive picture in mind, we consider the scenario in which the FE order is dominantly generated by a lattice distortion. For simplicity we consider a cubic ionic crystal with two atoms in the unit cell (the rocksalt structure of the IV-VI semiconductors, see Fig. 1). We label the two sublattices by b and r corresponding to the "blue" and "red" ions, which have equal in magnitude and opposite sign charges. Each sublattice has a corresponding phonon displacement field u r and u b . As usual, there are two modes: a gapless acoustic mode given by the sum u ac = (u r + u b )/2 and a gapped optical branch given by the difference u = u r − u b . Near the FE transition, the optical branch becomes nearly gapless and is described by the effective action Here c L and c T are the longitudinal and transverse phonon velocities, respectively, ω T is the phonon mass, which is the tuning parameter to the transition, and V is the lowest order symmetry allowed anharmonic correction to the phonon energy (where we have neglected additional anisotropic terms allowed by the cubic symmetry [37]).
The Coulomb term -The third quadratic term is that of the Coulomb potential: where ε ∞ is the bare dielectric constant, which accounts for the contribution of core electrons. This contribution is due to the transitions between the high-energy atomic configurations, and does not include the contributions from the lattice dynamics or electronic interband transitions close to the Dirac point.

B. Coupling terms
We now consider the coupling terms between the fields introduced in Eqs. (1)-(3). Electron-Coulomb coupling -We start with the coupling between the Dirac electrons and the Coulomb potential where ρ e = nψ n γ 0 ψ n is the electronic density.
Phonon-Coulomb coupling -The coupling of the ferroelectric phonon modes to the Coulomb potential follows from Eq. (4) by noting that the deviations of the "red" and "blue" ionic density from equilibrium are given by ρ r = ρ 0 − ∇u r and ρ b = ρ 0 − ∇u b . Given that the ionic charges are of equal magnitude and opposite signs, the coupling of the lattice to the Coulomb field φ is given by where Q is the ionic charge on "blue" sites (charge on "red" sites equals −Q).
Notice that the form of the coupling (5) implies that only the longitudinal phonon mode couples to the Coulomb field. We also point out that after integrating Eq. (5) by parts one gets a dot product between the polarization density P = Qu and the electric field E = −i∇φ. Therefore, this equation can also be viewed as the action of a dipole moment density in an electric field. Finally, in the case of a non-polar covalent crystal (e.g. elemental bismuth) all atoms in the unit cell are neutral, leading to a vanishing coupling Q = 0. Electron-Phonon coupling -We now consider the coupling between the Dirac electrons and the phonon modes. We write down this coupling from general symmetry arguments. The phonon mode u is a time-reversal invariant vector. Inspecting all possible local Dirac bilinears specified in Table I, we find that the only Dirac bilinear that forms a time-reversal symmetric vector and, thus, is allowed to couple to the phonon displacement field is ψγ 0 γ j ψ. Therefore, the corresponding coupling is given by In Eq. (1), we have assumed the Dirac cones occur at the inversion symmetric points in the Brioullin zone. In the case they do not, the coupling can also include interflavor scattering. We emphasize that, for simplicity, we consider a rotationally symmetric model in the main text. We discuss the possible effects of the cubic anisotropy in Appendix D.

III. RENORMALIZATION GROUP ANALYSIS NEAR THE CRITICAL POINT
We now use RG to analyze the theory introduced in the previous section. Summing up Eqs. (1)-(6) we have The coefficients Z ψ and Z u account for the renormalization of the dynamical terms.
Here, Λ 0 is a momentum UV cutoff corresponding to the scale at which electron dispersion can be considered linear, and l is "RG time". As the second step, we further rescale momenta, frequencies, and the long-wavelength parts of the fields according to to restore the UV cutoff e −δl Λ 0 back to Λ 0 . Above, z is the dynamical exponent, and η ψ , η u , η φ are engineering field dimensions. This rescaling leads to the tree-level RG flows of the couplings after coarse-graining by the factor e l (the argument l is suppressed for brevity): It should be mentioned that the choice of dynamical and field exponents is somewhat arbitrary here since it does not affect the flow of dimensionless coupling constants [44,45]. The special choice ω 2 T = 0, Q = e = λ = 0, η ψ = 5/2, η u = η φ = 3, and z = 1 makes the theory scale invariant. This is a non-interacting fixed point. At this fixed point, ω 2 T and Q are relevant perturbations, while e and λ are marginal at the tree level. Since ω 2 T is the tuning parameter for the FE transition, we will assume it small close to the critical point. In what follows, we focus on two distinct cases: The case of ionic crystals with Q = 0 and the case of covalent crystals with Q = 0.

Fixed point theory
Near the non-interacting fixed point introduced above, Q is relevant and, at the tree level, obeys the following RG equation: Thus, in the case of ionic crystals, Q grows rapidly to strong coupling. Therefore, we should first derive the effective low-energy theory with large coupling Q (of the order of UV cutoff) and then proceed to the RG analysis of the resulting theory. We can integrate out the Coulomb field φ, which generates the following terms The first term is the standard Coulomb interaction. The second term can be viewed as a phonon mass generated in the longitudinal sector (note that it is independent of the magnitude of the momentum). This mass generation is the well-known LO-TO splitting in ionic crystals [46]. Finally, the last term is the Fröhlich coupling between the longitudinal mode and electronic density.
The generated mass term for the longitudinal phonon mode is of the form ω L ≡ 4πQ 2 /ε ∞ + ω 2 T . Consider the limit of large Q, such that ω L /c L ∼ Λ 0 , where Λ 0 is the UV cutoff. We can further integrate out the massive longitudinal phonon mode and only keep the leading order terms in the 1/Q expansion. This procedure generates the standard dynamically screened Coulomb interaction where is the dynamical dielectric constant, which manifestly satisfies the Lyddane-Sachs-Teller relation [46]. Close to the critical point ω T → 0, which implies that the dielectric constant scales as ε(ω, q) ≈ ε ∞ ω 2 L /(ω 2 + c 2 L q 2 ) and diverges at low energies and momenta. Thus, the effective fine-structure constant α = e 2 /v F ε ∞ , which signifies the strength of the Coulomb interaction, becomes highly irrelevant and flows quickly to zero. It means that the Coulomb interaction is effectively screened by the longitudinal phonon mode. In the end, the FE critical point is controlled by the following effective field theory where is the projector to the plane transverse to q. The couplings V and λ here are weakly renormalized after integrating out the longitudinal mode. Before continuing we would like to make an important remark. We have arrived at a low-energy effective theory describing soft fluctuations of the phonon modes, where the only coupling to the charged fermions is the direct coupling Eq. (6). This is a somewhat counter intuitive result, as the phonon distortions generate huge dipolar moment with strong Coulombic forces. One might expect that these Coulomb forces are still dominant over the direct coupling we consider. We note, however, that the same Coulomb forces that cause the phonon to couple strongly to electrons also cause strong coupling of the phonons to themselves -the coupling that caused the generation of the mass term ω L and the LO-TO splitting. Thus, the polarization of the soft modes is precisely the part of the phonon that does not generate a dipolar moment. Thus, the only relevant coupling to the dynamical soft modes is Eq. (6).

One-loop RG analysis
Now we analyze the effective field theory for ionic crystals (14) within the one-loop RG approach. To get rid of the exponents η ψ , η u , z, which in principle can be chosen arbitrary, we focus on the dimensionless quantities which are independent of these engineering dimensions [44,45]. First, we derive coupled RG equations for the ratio of the phonon to electron velocities ν T ≡ c T Z ψ /v F Z u and the dimensionless electron-phonon coupling constant β ≡ λ 2 /4πc 2 T v F Z ψ (the details of the calculation can be found in Appendices A and B): The most important result that can be extracted from these equations is that the electron-phonon coupling β flows to the strong-coupling regime, see Fig. 2. Consequently, we conclude that the 3+1D ferroelectric quantum critical point in a Dirac semimetal considered in this paper is generically a strongly-coupled problem, even if the original UV value of the coupling constant is small. This may be contrasted with standard QED in three dimensions, where the flow of the coupling is towards weak coupling and the low-energy effective theory is the ?ree Dirac dispersion with renormalized parameters [47]. In the next section we discuss the possible superconducting instabilities resulting from this flow to strong coupling.
Our RG equations were derived under the assumptions of the zero Dirac fermion mass and Fermi energy, while the one-loop approximation is valid provided the coupling remains small. Given the flow to the strong coupling, it is important to understand what stops the RG flows. Here we estimate the scale at which β becomes of order 1 and defer the discussion of a finite Dirac mass/Fermi energy to Sec. IV. In realistic materials, the Fermi velocity is much bigger than the phonon velocity, thus, one can set ν T ≈ 0 in Eq. (16). Then, the equation for the flow of β can be readily integrated. Completely neglecting the mass of the soft mode, ω T ≈ 0, we find that β grows to ∼ O(1) at the RG scale l β = 3π/N β 0 , which corresponds to the momentum scale Here β 0 1 is the initial UV value of the coupling constant at the scale Λ 0 .
Another natural scale that serves as a cutoff for our RG equations is set by the flow of the (dimensionless) mass of the transverse phonon mode r ≡ ω 2 T /c 2 T Λ 2 0 , which determines the critical region: where γ ≡ V /c 3 T Z u is the dimensionless phonon-phonon interaction. Assuming that β and γ are small compared ) and the electron-phonon coupling β0. The red dashed line separates the two regions Λr > Λ β and Λr < Λ β , corresponding to weak-and strong-coupling regimes, respectively. In the former region, the finite mass of phonons cuts off the RG flow before β reaches strong coupling, and the theory with renormalized parameters remains weak-coupled. In the latter case, β flows to strong coupling before the system leaves the critical region. This regime is characterized by strong electron-phonon coupling and requires further study. The insets in each region schematically depict the dispersion close to the four L points. Every Dirac cone in the paraelectric phase splits into two Weyl points in the ferroelectric phase.
to the UV value r 0 , the solution of this equation with the exponential accuracy reads as r ∼ r 0 e 2l . The critical regime is determined by the condition r 1, which corresponds to the RG scale l r ≈ (1/2) ln(1/|r 0 |), or, equivalently, momentum scale [48] If Λ β > Λ r , the theory flows to the strong coupling regime before the phonon mode gets massive. Our RG equations are only applicable then down to Λ β . In the opposite case, Λ r > Λ β , the RG flow should be stopped at Λ r , where the transverse phonon mode becomes massive and can be integrated out. At this scale, the system leaves the critical regime, while the coupling between phonons and fermions still remains weak. The corresponding phase diagram is shown in Fig. 3. We will consider the latter case in more detail in the next section in context of superconductivity.
Finally, we discuss the flow of the dimensionless phonon-phonon interaction γ ≡ V /c 3 T Z u , which corresponds to the anharmonicity of the lattice oscillations: This equation, again, can be easily analyzed in the physical case ν T ≈ 0. Then, since β is a marginally relevant parameter, γ eventually also flows to strong coupling. It is straightforward to show, however, that this flow does not introduce any new cutoff, as γ can reach order 1 no sooner than at Λ β given by Eq. (17), which is realized in the large-N limit (i.e., when the term proportional to ∝ γ 2 on the right-hand side of Eq. (20) can be neglected). It is also interesting to note that sufficiently large ν T in Eq. (20) can, in principle, drive γ negative, thus indicating a first-order transition into the ferroelectric state.
Since we consider ν T ∼ 1 hardly realizable in real physical systems, we do not study this possibility in detail here.
Another interesting result that can be inferred from the RG equations is the flow of the electron and phonon velocities (here we fix the dynamical critical exponent z = 1): We see that one of the physical properties of the ferroelectric critical point in Dirac materials is the reduction of the velocities under RG for both the transverse phonon modes and the Dirac fermions. Furthermore, as is shown in Fig. 2, for N = 1, the velocity ratio ν T flows to one of two possible values ν 0 = 0 or ν 1 = 1, depending on whether the initial value of ν T is smaller or larger than ν m = t/3 2/3 − 2/(3 1/3 t) − 1 ≈ 0.18, respectively, with t = (18 + 2 √ 87) 1/3 . If N > 1, the flow is always towards ν 0 = 0.
So far we only considered a rotationally symmetric model with isotropic electron and phonon velocities. For N > 1, however, there is no symmetry that forbids anisotropic terms that manifest the symmetry of the underlying lattice. Nevertheless, the accounting for these terms does not modify main qualitative results described above. Hence, we focus on the isotropic case for the rest of the paper for simplicity, and defer the discussion of possible anisotropies to Appendix D.

B. Covalent crystals (Q = 0)
Now we perform similar RG analysis for covalent crystals exemplified by elemental bismuth. While the main qualitative results, such as flow to strong coupling, in this case are the same as for ionic crystals, certain important differences should be discussed. In particular, as mentioned above, the optical phonon distortion u generates a negligible amount of polarization in covalent crystals. Consequently, the effective theory for these materials is described by Eq. (7) with Q = 0. As a result of this important difference, the argumentation of Section III A about the screening of Coulomb interaction by longitudinal phonons no longer holds. Instead, one should keep track of the flows of the parameters ε ∞ and e, in addition to those considered in Eq. (16). Focusing again on dimensionless parameters that do not depend on engineering dimensions η ψ , η u , η φ , and z, we find the following system of coupled one-loop RG equations: Since both the longitudinal and the transverse phonon modes become massless at the transition in covalent crystals, they should be treated on equal footing. Consequently, one could in principle consider two (not independent) dimensionless couplings β = λ 2 /4πv F c 2 T Z ψ and β = λ 2 /4πv F c 2 L Z ψ , which quantify the electron-electron interaction strength mediated by the transverse phonons and the longitudinal phonons, respectively. It is straightforward to show, however, that, in the physical limit ν T ∼ ν L 1,β is marginally irrelevant, while β flows to strong coupling. Indeed, in this limit, first two equations of (22) take form while analogous equation forβ would read as Similarly to the case of ionic crystals, electron-phonon coupling β flows to the strong-coupling regime, while the Coulomb interaction α becomes suppressed under RG. The important difference, however, is that now α is only marginally irrelevant and flows to zero much slower. The reason for this difference is that the Coulomb screening in covalent crystals is due to interband (between particle and hole bands) transitions, which is much weaker than the screening due to lattice polarization in ionic crystals. The flow of the phonon velocities can also be easily studied in the limit ν T ∼ ν L 1. Analogously to ionic crystals, ν T flows to zero in this regime. The flow of ν L , on the other hand, is sensitive to the number of flavors N and to the initial conditions, as well as to the scale that stops RG. For instance, at sufficiently large N , ν L is increased under RG.
Finally, the flows of the phonon-phonon coupling γ and the phonon mass r are qualitatively similar to the case of ionic crystals, so we do not consider them here in detail.

IV. SUPERCONDUCTIVITY
In the previous section we have analyzed the RG flow of the electron-phonon coupling near a ferroelectric quantum critical point. We found that, generically, the critical point is unstable and flows to strong electron-phonon coupling, while the Coulomb interaction flows to weak coupling. As a result, we anticipate that the effective electron-electron attraction mediated by the ferroelectric phonon modes will become dominant over the Coulomb repulsion. Hence, the natural next step in our study is to apply this result to superconductivity.
We emphasize that in our scenario, both for ionic and covalent crystals, the enhancement of the attractive interaction over the Coulomb repulsion does not require finite electron density, in contrast to the Anderson-Morel theory. Nonetheless, this does not imply that the superconducting transition temperature does not depend on density. At least at weak coupling, β 1, a finite density of states is essential to obtain a finite T c . Therefore, we will now relax our previous assumption about the Fermi energy exactly at the Dirac point and assume a finite Fermi momentum k F . As before, we separately consider the cases of ionic and covalent crystals. We also focus on the paraelectric side of the transition, i.e., consider systems possessing both time-reversal and inversion symmetries in the normal state.

A. Ionic crystals
As we have shown in Sec. III, one can define two scales Λ β and Λ r given by Eqs. (17) and (19), which denote the divergence scale of the electron-phonon coupling β and the phonon mass r, respectively. When Λ r > Λ β , r diverges first and thus the flow is terminated before β reaches strong coupling (this regime is denoted by the shaded regions in Fig. 3). In what follows we consider this weak coupling limit, where the BCS approach is applicable, and leave the strong coupling regime Λ β > Λ r for a future work.
The additional scale we have introduced, k F , can, in principle, also put the flow to a halt when the running scale Λ(l) = Λ 0 exp(−l) becomes of order k F . Thus, depending on the ratio between k F and Λ r , one may again consider two cases. The first case, k F Λ r , is close to the standard Anderson-Morel scenario with the phononassociated scale (ω T /c T )| l=0 being smaller than k F and we do not consider it here in detail. Since we are interested in understanding superconductivity at very low density, we focus on the opposite case Λ r k F . In this limit, the screening of the Coulomb repulsion by longitudinal phonons occurs well above the Fermi scale, as disscused below Eq. (13), and we obtain a Fermi liquid with static phonon-mediated attraction.
To obtain an effective low-energy interaction, we allow the system to flow according to the RG equations derived in Sec. III A until it reaches the scale Λ r . We then use Eq. (14) to integrate out the transverse phonon mode, which is massive at this scale, with the effective propagator that can be considered frequency-and momentumindependent. This procedure results in the attractive interaction Hamiltonian where the effective interaction constant is obtained from Eq. (16) in the limit of ν T 1. To make the analysis similar to the conventional BCS at this point, we write Eq. (25) in the Hamiltonian formalism (and use ψ † instead ofψ). This became possible since at the scale Λ r the phonon-mediated interaction can be considered static, ω ω T , analogously to BCS theory.

Projection onto the conduction band
Now we analyze the superconducting instabilities due to interaction (25). We assume that the Fermi energy ε F = v F k F is much larger than the superconducting gap, ε F ∆, hence, the conventional weak-coupling BCS-like treatment is applicable. In this case, it is convenient to project all operators onto the conduction band, thus significantly simplifying the model by reducing it from the original four-orbital to effective two-orbital. In the paraelectric phase, the only case we consider in this Section, both time-reversal and inversion symmetry are present in the normal state, hence, all energy bands remain double degenerate even in presence of strong spin-orbit coupling. The electron states are characterized by two-component T . In the presence of spin-orbit coupling, however, components c 1,2 are not spin eigenstates anymore, but rather eigenstates in some band basis. The choice of this basis, however, is not unique. For concreteness, we choose the so-called manifestly covariant Bloch basis (MCBB), in which [c 1 (k), c 2 (k)] T transforms as an ordinary SU (2) spin-1/2 [49]. To find this basis, we diagonalize Hamiltonian which corresponds to Eq. (1), and choose the band eigenstates to be fully spin-polarized along the z-axis at the origin of the point group symmetry operations (see also Refs. [50] and [51] for more details). The eigenvectors a 1 (k) and a 2 (k) in the MCBB that correspond to the states in the conduction band are given by where η = ±1 corresponds to the electron/hole band, respectivly, and we definedk ± = (k x ± ik y )/k. The mapping onto the MCBB then simply implies the transformation ψ(k) → a 1 (k)c 1 (k) + a 2 (k)c 2 (k), and can schematically be written as ψ(k) = Q η (k)c(k), where Q η (k) is a projector onto MCBB. It is straightforward to show then that the Dirac bilinear γ j , which couples to a soft ferroelectric mode, projects onto where lmj is the Levi-Civita tensor, and σ j here are Pauli matrices acting in MCBB.

Pairing channels and transition temperature
To demonstrate the superconducting instabilities, we now decompose interaction (29) into pairing channels, analogously to how it has been done in Ref. [50]. The time-reversal invariant superconducting order parameter generally takes form where, again, ε βγ is the Levi-Civita symbol. In systems with strong spin-orbit coupling, spin S and angular momentum L are not good quantum numbers. Instead, in systems with O(3) symmetry considered here, all possible orders are characterized by the total angular momentum J = L+S. As was shown in Refs. [49] and [50], the formfactors F αβ (k) up to order J = 1 have the form shown in Table II. L = 0 state F 0 corresponds to the conventional s-wave pairing with J = 0, while L = 1 sates are odd-parity p-wave, and transform as a pseudoscalar (F 1 with J = 0) and a vector (F j 2 with J = 1) under the symmetry operations.
Next, we restrict the effective interaction (29) to the Cooper channel with the zero total momentum by keeping terms with k = −k only. Focusing on the states near the Fermi surface, |k| ≈ |k | ≈ |k + q| ≈ |k − q| ≈ k F , it is straightforward to decompose Eq. (29) into the pairing channels F n [50]: where coefficients a T n are listed in Table II. The ellipsis on the right-hand side of Eq. (31) denote terms with J > 1 [52]. The contribution from these terms is numerically small, and we do not consider it in this paper.
Up to order J = 1, only two channels are attractive and lead to a superconducting instability: the scalarF 0 with a T 0 = 1 and the vectorF j 2 with a T 2 = 1/2 [53]. We thus, conclude that pairing in the s-wave channel is the most dominant superconducting instability. However, the projector significantly reduces T c in this channel without modifying it in the p-wave channel with n = 2. As we demonstrate in the next section when considering the covalent crystals, this difference may play significant role leading to an odd-parity superconductor if the Coulomb interaction is not vanishingly small.
The transition temperature T c can be estimated from Eq. (31) using the usual gap equation [54]: In case of the most attractive s-wave channel, it takes form where ν 0 = k 2 F (l r )/2π 2 v F (l r ) is the density of states at the Fermi energy per one spin projection per one Dirac node, and all quantities entering it are taken at the RG scale l r . We emphasize that the upper cutoff in Eq. (32) is not the phonon frequency, as in the standard BCS theory, but given by the Fermi energy. This situation is somewhat analogous to the superfluidity in a chargeneutral Fermi liquid, studied in Ref. [55]. We estimate transition temperature from Eq. (32) as The parameters k F , c T 0 , and ω T 0 are the original (UV) values of the Fermi momentum, phonon velocity, and phonon mass, respectively, while β * is renormalized according to Eq. (26), and we used k F (l) = k F e l . Hence, the proximity to the ferroelectric critical point leads to a significant enhancement of T c . To emphasize this point, we rewrite Eq. (34) in the form where T c0 ∼ ε F exp −πω 2 T 0 /2k 2 F c 2 T 0 β 0 is the estimate for a transition temperature that we would obtain without taking into account the critical nature of the ferroelectric fluctuations. We see that, even within the weak-coupling approximation β 0 ln(Λ 0 /Λ r ) 1, we obtain huge enhancement of the transition temperature by a factor (1/r 0 ) λ . This result is to some extent similar to the enhancement of T c by the critical nematic fluctuations obtained in Ref. [30].
Finally, we estimate the temperature that would correspond to a transition into the p-wave superconducting state T cp : (36) An additional factor 3/2 in the exponent appears due to the averaging over the directions of vector k in Eq. (32). T cp is exponentially smaller than T c , and, consequently, p-wave superconducting phase seems unreachable within the present scenario. However, we demonstrate in the next section that the presence of the repulsive Coulomb interaction can, under certain conditions, suppress swave channel and drive a system into the odd-parity superconducting state.

B. Covalent crystals
Our analysis of superconductivity in covalent crystals is similar to the ionic case. There is, however, two im-portant differences. First, the longitudinal phonon mode also becomes soft at a ferroelectric transition, consequently, there will be an additional contribution to the effective electron-electron interaction mediated by a longitudinal mode. Second, the Coulomb repulsion is not screened by the lattice polarization, but only by the interband transitions. As we showed in Sec. III B, the correspondent coupling constantsβ and α are marginally irrelevant, see Eqs. (23) and (24). They flow to zero only logarithmically upon RG and, thus, should be taken into account in the weak-coupling regime we are considering here. As we show below in Sec. IV B 1, the inclusion of the Coulomb interaction allows, upon proper tuning of the coupling constants, to realize a p-wave superconductor.
The effective electron-electron interaction due to longitudinal phonons projected onto the conduction band has the same form as Eq. (25), but with the substitution P ij (q) → δ ij − P ij (q) and β * →β * , whereβ * is given by (see Eq. (24)) The decomposition into the pairing channels has the form (31), with β * substituted byβ * and coefficient a L n are listed in Table II. We see that the interaction mediated by the longitudinal phonons also favors s-wave pairing, hence, its only effect is to modify the expression for T c accordingly. The inclusion of the Coulomb interaction, on the other hand, may have more dramatic consequences, leading, under certain conditions, to the p-wave superconductivity in covalent crystals.

Possibility of p-wave pairing
To demonstrate how the Coulomb interaction may result in the p-wave superconductivity, we generalize our analysis for the case of a finite Dirac mass m in Eq. (1). Again, we focus on the regime with m/v F , k F Λ r , so the RG flow is not affected by the finite mass/Fermi energy, and stops at the same scale Λ r , while the ratio m/v F k F can be arbitrary.
In case of a finite mass, the eigenvectors in the MCBB (27) are generalized as where we only consider states at the Fermi surface, |k| = k F , and defined β ± = (1/2) 1 ± m/ε F with (40) also assumes the Fermi energy inside the electron band, while the expression for the opposite case is obtained by the substitution The effect of a finite mass on the phonon-mediated part of the interaction is rather simple: all coefficient a n in Eq. (31) (for both longitudinal and transverse phonons) should be replaced by a n → (v 2 F k 2 F /ε 2 F )a n , and the density of states ν 0 in Eqs. (34) and (36) should be modified according to ν 0 = ε F k F /2π 2 v 2 F . Again, all the parameters entering ν 0 (as well as mass m) are taken at the scale Λ r , and we suppress the argument l r for brevity.
The decomposition of the Coulomb repulsion is more subtle.
Because of its long-range nature, the qdependence of the interaction must also be taken into account. Taking the simple Thomas-Fermi approximation and projecting onto the MCBB, we find: where q 2 T F = 8πN α * ν 0 v F is square of the Thomas-Fermi wavevector, and α * is given by Again, all quantities entering q T F are taken at the RG scale l r . Focusing again on the states at the Fermi surface only, we decompose the Coulomb interaction (41) into the pairing channels: where the ellipsis stands for the terms with J > 1.
The expression for functions f n (x) are presented in Appendix E. Summing up contributions from the transverse and longitudinal phonon modes, Eqs. (31) and (39), and direct Coulomb repulsion (43), the decomposition of the total effective electron-electron interaction into the pairing channels has form where coefficients a T n and a L n are presented in Table II, and functions f n (x) are listed in Appendix E. Equation (44) is a direct generalization of Eq. (31) for the case of covalent crystals and finite Dirac mass. The expression for T c in the attractive pairing channels should also be modified accordingly.
In general, functions f n (x) have rather complicated form. However, to demo?strate how the p-wave superconductivity may appear, it is sufficient to consider the limit of a very low density, v F k F ε F . Assuming further that α * is not too small, we find It is clear from the above expression that, as long as α * ε F v F k F , the s−wave channel is much more suppressed by the Coulomb repulsion than the p−wave channel, f 0 f 1,2 . We further assume that the coupling constants β * andβ * are renormalized significantly enough, such that β * ∼ 1 andβ * ≈ 0. It means that the system is on the verge of entering the strong-coupling regime, while the contribution from the interaction mediated by the longitudinal phonons can be neglected. Then, the ratio between the phonon-mediated attraction and the Coulomb repulsion in the s-wave channel can be rudely estimated as We see that, because of the small factor v F k F /ε F 1, the Coulomb repulsion significantly exceed the attraction due to phonons, thus completely suppressing superconductivity in this channel.
On the other hand, the analogous estimate for the vector-type p-wave pairing channel F i 2 , which, according to Table II, is also attractive if only the transverse phonons are considered, gives where we explicitly restored the argument l r . We see that, unlike the s-wave, the above expression does not have the smallness v F k F /ε F . Consequently, assuming that α * , β * ∼ 1, the ratio (47) can be of order 1 provided the smallness k F (l r )/Λ 0 = k F /Λ r ∼ k F c T 0 /ω T 0 is compansated by a large numerical prefactor and a large number of Dirac cones N .
While the prerequisites for the p-wave superconductivity in the described above scenario impose a lot of constraints on the parameters entering the problem, we believe that the qualitative conclusion about the possibility of the p−wave pairing from the interplay between the electron-phonon attraction and the Coulomb repulsion remains valid. A careful quantitative study of this problem, however, most likely requires methods suitable for analyzing the strong-coupling regimes, such as the numerical Quantum Monte-Carlo or the analytical Migdal-Eliashberg approaches. We leave these problems for the future study.

V. CONCLUSIONS
We have studied the ferroelectric quantum critical point in low-density Dirac materials for the cases of ionic and covalent crystals. We derived a general low-energy effective field theory that includes the interaction between soft phonon modes and electrons, as well as direct Coulomb repulsion. We showed that, in case of polar crystals, the dominant interaction between electrons is mediated by the transverse phonon mode, while the Coulomb repulsion is screened by the lattice. Using RG analysis, we demonstrated that the effective low-energy theory flows to a regime with strong electron-phonon coupling. Similar results were obtained for covalent crystals. In the latter case, however, the screening of the Coulomb repulsion is much weaker due to the lack of lattice polarization. We further demonstrated that the proximity to the FE critical point significantly enhances superconductivity. Finally, we showed how the interplay between phonon-mediated attraction and Coulomb repulsion in covalent crystals can, in principle, lead to p-wave superconductivity.
It is interesting that the problem of a Dirac semimetal undergoing a ferroelectric transition generically flows to strong coupling even in the absence of a finite Fermi surface. This problem can be addressed using determinental quantum Monte-Carlo, since it does not suffer from a sign problem. It will be informative to study if the superconducting transition survives in the strong-coupling regime, or it is destroyed by strong critical fluctuations.
We also expect that the strong-coupling regime close to the critical point will have experimental consequences. For example, the high-accuracy measurements of the phonon dispersion may reveal a 2k F Kohn anomaly. Furthermore, strong coupling will significantly enhance electrons' scattering, leading to a strong temperature dependant resistivity.
Considering our results in a broader context, we expect this mechanism to be relevant to all low-density superconductors that posses a nearly touching conduction and valence bands. As explained in the introduction, this applies to almost all low-density superconductors, including bismuth, YPtBi, PbTe, SnTe, Sr x Bi 2 Se 3 , Ge and Sr 3−x SnO. We also find it interesting to speculate that this mechanism may be relevant to the high-T c lowdensity superconductor FeSe, which possesses a Diraclike dispersion close to the M -points in the Brillouin zone.

VI. ACKNOWLEDGMENTS
We thank Patrick Lee, Liang Fu, Meng Cheng and Kamran Behnia for illuminating discussions. ZB acknowledges support from the Pappalardo fellowship at MIT.
Here we list the one-loop diagrams that contribute to our renormalization group equations. We start from the action for the covalent crystals by including both the transverse and longitudinal optical phonon modes. The result for the polar case can be deduced by restricting to the diagrams with only transverse phonons (formally just take e → 0, c L → ∞ limit).
In this paper, we use the standard momentum-shell RG scheme by the integrating out an infinitesimal momentum shell Λe −δl < q < Λ, and all frequency −∞ < ω < ∞. Fermion self energy -Both the optical phonon and the coulomb field contribute to the fermion self energy. The contributions from the phonon fields come from both the transverse phonon field and the longitudinal phonon field.
The self energy from the Coulomb interaction reads Optical phonon self energy -The polarization of the phonon field given by the fermion bubble diagram is: Interestingly, we notice that the fermion bubble diagram renormalizes the velocity of the transverse mode and the longitudinal mode in opposite ways. In this sense, the fermions can make the two modes very different. When considering the ionic crystal case, we simple omit the renormalization to the longitudinal mode since it has a large gap already.
Additionally, the self interaction of the boson field also generates a self energy correction, which is given by Coulomb field self energy - Electron-phonon vertex correction -Each vertex correction has two contributions, one where the boson exchanged in the loop is the same boson of the vertex and one where it is the other bosonic field. For example the vertex correction to the coupling of the phonon is given by Coulomb vertex correction -We have corrections to the Coulomb vertex: Phonon interaction vertex correction-Finally, the four phonon vertex is renormalized by the phonon bubble diagrams and also the fermion box diagram. Over all we find the correction is given by As explained in Section III A, in the ionic case the Coulomb field and longitudinal phonons develop masses similar to the Higgs mechanism and as a result the Coulomb interaction and longitudinal phonons becomes irrelevant for our renormalization group study. Thus, the RG equations for the ionic case derive from the calculations in Appendix A by setting e = 0 and c L → ∞. This leads to the following RG equations The RG equations for the covalent case can be readily obtained from Eqs. (A2)-(A11): Appendix D: Comments about crystal anisotropy In the analysis presented in the main text we have considered a fully rotational invariant system. In a realistic crystal, however, there are always anisotropies. In what follows we will discuss such anisotropies in ionic crystals with cubic symmetry.
A cubic anisotropy has two important effects that are relevant to the flow of β. First, the dsipersion of the soft modes (2) will include the anisotropy term S a u = −κ d 4 q q 2 j u 2 j We negelct this term in what follows, i.e. we assume κ = 0. We also mention a recent comment where the effect of this term on the polarization of the ferroelectric modes was computed perturbatively [38]. The second effect of crystal anisotropy, which we will consider, appeares in the electronic dispersion, which is relevant only for N > 1. Let us consider N = 4 where the four Dirac points occur on the boundary of the BZ at the L-points. In this case the Dirac dispersion term Eq. (1) assumes the form [39] where z-direction is defined differently for each Dirac point. It corresponds to the line connecting the Γ-point to each of the L-points. Additionally, the coupling of each of these Dirac points to the phonons Eq. (6) is also becomes modified where λ z is the coupling to a distortion along the line and λ is the coupling transverse to it. The RG procedure described in Sec. III can still be performed analytically, although the expressions become rather lengthy. We also have two additional dimensionless parameters v z /v F and λ z /λ that flow under RG.
The crucial point is that the main result of this paper, the flow of β towards strong coupling remains unchanged. To demonstrate this we plot the β-function for β in the anisotropic case normalized by the isotropic case Eq. (16) for ν T = 0.1 in Fig. 5. Here the ratio is plotted for three different values of λ z /λ as a fucntion of v z /v F . We find that in all cases the β-function is positive. Also note that the regime of v z < v F and v z > v F should not be considered on equal footing, since in the latter case the density of states at the Dirac point is enhanced while in the opposite limit it is decreased. Also, note that we did not perform a detailed study of the multidimensional flow in the four dimensional space of all four parameters.
The asymptotic behavior of these expression in the limit k F v F ε F = m 2 + k 2 F v 2 F is presented in the main text.