Quasiparticles near domain walls in hexagonal superconductors

Quasiparticles near domain walls in hexagonal superconductors

S. P. Mukherjee and K. V. Samokhin Department of Physics, Brock University, St. Catharines, Ontario, Canada L2S 3A1
Abstract

We calculate the energy spectrum of quasiparticles trapped by a domain wall separating different time reversal symmetry-breaking ground states in a hexagonal superconductor, such as UPt. The bound state energy is found to be strongly dependent on the gap symmetry, the domain wall orientation, the quasiparticle’s direction of semiclassical propagation, and the phase difference between the domains. We calculate the corresponding density of states and show how one can use its prominent features, in particular, the zero-energy singularity, to distinguish between different pairing symmetries.

pacs:
74.20.-z, 74.55.+v

I Introduction

The presence of domain walls (DWs) in a superconductor is a direct evidence of an unconventional nature of the pairing, because the DWs can only appear if there are two or more distinct degenerate ground states, which transform one into another by some discrete symmetry operations, e.g., by time reversal. This is possible if the superconducting order parameter has more than one component, i.e., corresponds to either a multidimensional irreducible representation (IREP) of the crystal point group or to a mixture of different one-dimensional (1D) representations. The chiral -wave state, which is realized, for example, in SrRuO (Ref. Mackenzie, ), is a well-known example of a system in which DWs are believed to play a prominent role. Strong evidence of the superconducting states with broken time reversal symmetry (TRS) has also been reported in URuSi (Ref. Schemm, ), UPt (Ref. UPt3-review, ), SrPtAs (Refs. SrPtAs-exp, and Fischer, ), and PrOsSb (Refs. PrOsSb-exp, and AC07, ). Various TRS-breaking states have been proposed theoretically in BaKFeAs (Ref. spis, ), doped graphene (Ref. NLC, ), undoped bilayer silicene (Ref. Liu2, ), NaCoOHO (Ref. KPHT, ), and other materials. Superconducting DWs can be created in these systems, e.g., due to the nucleation of the order parameters of opposite chirality in different parts of an inhomogeneous sample.

It is well known that a superconducting DW can trap quasiparticles in its vicinity, creating the Andreev bound states (ABS), see, e.g., Ref. p-wave-ABS, . The energy of these states is inside the bulk gap and their very existence can be explained by topological arguments, see Refs. Volovik, and Bernevig, . The ABS contribution to the tunneling density of states (DOS) can be easily separated from that of the bulk quasiparticles and can, therefore, be used to prove the DW presence. Moreover, the ABS spectrum is sensitive to the gap structure in the bulk of the domains, which allows one to confirm or rule out certain pairing symmetries.

We focus on the case of a three-dimensional (3D) hexagonal superconductor with the crystallographic point group , which describes UPt. The quasi-two-dimensional tetragonal case, which is applicable to SrRuO and the iron-based superconductors, was previously studied in Ref. Samokhin1, . The heavy-fermion superconductor UPt has a complicated phase diagram, with two distinct phases (called and phases) even in the absence of external magnetic field, see Ref. UPt3-review, for a review. A variety of thermodynamic and transport measurements have revealed an unconventional superconducting state with nodal excitations. Although there is still no general consensus on the pairing symmetry in UPt, the most promising candidate model, which has recently received further support from the Josephson interferometryUPt3-phase () and the polar Kerr effectUPt3-Kerr () experiments, is based on the two-dimensional (2D) IREP of the point group . The corresponding order parameter is real in the high-temperature phase (at 0.45K 0.5K), and complex, i.e. TRS-breaking, in the low-temperature phase (at 0.45K).

Our goal is to study the quasiparticle tunneling features which are uniquely associated with the DWs and can be used to probe the symmetry of the superconducting order parameter, e.g., in the phase of UPt. We analyze the TRS-breaking states corresponding to all IREPs of that support the formation of DWs. The paper is organized as follows. In Sec. II, we derive a general expression for the ABS energy in the semiclassical (Andreev) approximation. In Sec. III, we calculate the ABS spectrum and the corresponding contribution to the DOS separately for each of the four possible TRS-breaking states. The summary of our results is presented in Sec. IV. Throughout the paper we use the units in which .

Ii Andreev bound states

We consider a hexagonal superconductor described by the point group , in zero magnetic field. The axis is along the sixfold symmetry axis and the plane coincides with the basal plane. The electron band dispersion is assumed to be , where is the effective mass, with generalization to a more general, e.g., ellipsoidal, case being straightforward. The superconductor is divided into two semi-infinite superconducting domains by a planar DW.

Since the scale of the order parameter variation in the DW is much greater than the inverse Fermi wavevector, we can use the Andreev approximation,And64 () in which the Bogoliubov quasiparticles propagate along the semiclassical trajectories characterized by the Fermi-surface wavevectors ( and are the spherical angles, with the polar axis directed along the positive axis). In the semiclassical approximation, the gap function is given by a spin matrix, which depends on the position and the wavevector : for singlet pairing, and for triplet pairing. We consider only 2D IREPs of , therefore the gap functions can be written in the following form (Ref. Mineev, ): and , where , are the order parameter components, and the basis functions satisfy , . The direction of the spin vector is assumed to be fixed along by the strong spin-orbit coupling of electrons with the crystal lattice.z-axis ()

In general, the quasiparticle wave function has four components, corresponding to the Nambu (electron-hole) and spin degrees of freedom, but in the models considered in this paper the spin channels are decoupled and the equations are reduced to a two-component form. For each spin projection, the wave function is a product of a rapidly oscillating plane wave and a slowly varying Andreev envelope function , which satisfies the equation

(1)

Here is the Fermi velocity and is the gap function sensed by the quasiparticles as they propagate along the semiclassical trajectory defined by : for singlet pairing, and for triplet pairing.

We use the sharp DW model, in which the gap function is described by two different complex constants in the two domains along each semiclassical trajectory:

(2)

where is normal to the DW plane. The angular dependence of the gap function is different for different IREPs of the point group, see Sec. III below. The solution of Eq. (2) which is exponentially localized near the DW has the form , where

From the continuity of the wave function across the DW, we obtain the following equation for the Andreev bound state (ABS) energy:

(3)

where . The solution of Eq. (3) is straightforward, see, e.g., Ref. Samokhin1, , and we find that for each direction of semiclassical propagation satisfying the condition

(4)

there is exactly one ABS, whose energy is given by

(5)

where and

(6)

In all cases studied in this work, we have , therefore and the condition (4) is satisfied for every direction. It is straightforward to show that the ABS energy (5) is inside the bulk gap, i.e., . The dependence of the ABS energy on the direction of semiclassical propagation is not continuous, showing abrupt changes when either or the Fermi velocity projection on the DW normal change their signs, see Appendix A.

The sharp DW model (2) can be justified by the following argument. The ABS wave function is exponentially localized on both sides of the DW, with the characteristic scales given by . The sharp DW approximation is legitimate for those directions of semiclassical propagation for which the DW width is smaller than . This condition is strongly angle-dependent and, in particular, fails for the trajectories corresponding to . However, for such trajectories the Andreev approximation itself is not applicable. For most directions of , one can use the following estimate: , where is a characteristic value of the gap and is the superconducting correlation length.

The quantity of interest is the DOS of the ABS’s, which can be measured in tunneling experiments. We consider two orientations of the DW, with the normal vector either parallel or perpendicular to the basal plane. In the first case, assuming , the order parameter depends only on , the momentum components parallel to the DW are conserved, and the DOS per unit DW area for both spin projections has the following form:

(7)

where and are the system’s dimensions in the plane. Derivation of this expression is outlined in Appendix B. Taking the thermodynamic limit and changing the integration variables to the spherical angles, we obtain:

(8)

where is the normal-state DOS in 3D at the Fermi surface per one spin projection. The ABS energy, see Eq. (5), has the following form:

(9)

In the case of , the order parameter depends only on , and we obtain the following expressions for the DOS per unit DW area for both spin projections:

(10)

and for the ABS energy:

(11)

Due to the electron-hole symmetry, , so that below we calculate the DOS only for .

Iii Andreev bound states spectra

The point group has twelve IREPs, six even and six odd, of which eight are 1D and four are 2D. The formation of DWs is possible only for those superconducting classes which are degenerate with respect to some discrete symmetry.VG85 (); Mineev () Since the 1D IREPs cannot support DWs, we focus on the 2D IREPs. We consider only the TRS-breaking chiral states, with the order parameters given by . The momentum dependence of the corresponding gap functions is listed in Table 1. Note that, due to the similarity of the basis functions, our results for the IREPs and are also applicable to a tetragonal superconductor with the point group .

IREP ,
Table 1: The momentum dependence of the singlet () and triplet () gap functions of the chiral states corresponding to the 2D IREPs of , for a strong spin-orbit coupling. The singlet gap functions correspond to the even IREPs and , while the triplet gap functions correspond to the odd IREPs and .

iii.1 representation

For the chiral -wave state corresponding to the IREP , we find from Table 1 the following expressions for the gap functions in the two domains: and . Here is the Josephson phase difference between the domains, which has to be included in order to satisfy the current conservation across the DW, see Ref. Sam12, and also below. Its value depends on the microscopic details of the system, but here we regard it just as an additional phenomenological parameter. In terms of the spherical angles the gap functions become

(12)

It follows from Eq. (6) that . Below we calculate the ABS energies and the DOS for both orientations of the DW.

. We obtain from Eq. (9) the following expression for the ABS energy:

(13)

It is shown, for and , in the upper panels of Fig. 1. For general , the energy is discontinuous at and also at , see Appendix A. The ABS energy has two lines of zeros at , which correspond to the quasiparticle trajectories with . These zeros in the ABS dispersion have a topological origin, see Sec. III.5 below. There are also two point zeros in the ABS energy at the poles of the Fermi surface, i.e. at and . However, for these directions (corresponding to the “grazing” trajectories, parallel to the DW) one has , and the Andreev calculation resulting in Eq. (5) is not applicable.

The quasiparticle DOS is given by Eq. (8) and can be found analytically for and . Since the calculation is similar in both cases, here we outline it only for , when we have

Since we consider only positive energies, the second delta function does not contribute to the integral and we obtain:

where . Evaluating the last integral,Prudnikov () we finally arrive at the following expression:

(14)

In a similar fashion, we obtain:

(15)

for . The DOS curves corresponding to Eqs. (14) and (15), normalized by , are shown in the bottom panels of Fig. 1. We can see that the overall magnitude of the ABS contribution to the DOS (per unit area) is of the same order as in the normal state, since .

. Using Eq. (11) we obtain for the ABS energy:

(16)

It is discontinuous at and also at , see Appendix A.

One can show that the current conservation requires that in the lowest order of the Ginzburg-Landau (GL) gradient expansion. The superconducting current can be obtained in the standard fashion from the gradient terms in the GL free energy density. Since the order parameter components depend only on , the gradient energy has the form , in the notations of Ref. Mineev, . Replacing the gradients by the covariant derivatives, , and varying with respect to the vector potential , we obtain for the superconducting current: . We use the constant-amplitude approximation for the order parameter:

where is the common (or Josephson) phase of the order parameter components and is the relative phase, satisfying . Then the supercurrent becomes . It follows from the current conservation that has a constant value, which is fixed by external sources. Setting , one obtains , and, therefore,

(17)

It is easy to see that this result holds for the chiral states corresponding to all four 2D IREPs. In contrast to the case of (Refs. Sam12, and Samokhin1, ), the Josephson phase difference between the domains for takes a universal value , i.e., does not depend on the coefficients in the GL expansion. This last conclusion can be invalidated by the inclusion of higher-order gradient terms and going beyond the constant-amplitude approximation. However, one can see from the way the angle enters Eq. (16) that the ABS dispersion for different is obtained from that for by simply translating the latter along the axis by . Therefore, the DOS, see Eq. (10), does not actually depend on . The calculation is similar to the case and the final result has the following form:

(18)

The ABS energy for and the DOS for any are shown in Fig. 2.

Figure 1: (Color online) The ABS energy as a function of the direction of semiclassical propagation (upper panels) and the corresponding DOS (lower panels), for and , in the case of the chiral -wave state () and . The grazing trajectories (for which ) are shown at and by horizontal dotted lines.
Figure 2: (Color online) The ABS energy as a function of the direction of semiclassical propagation for (upper panel) and the corresponding DOS (lower panel), in the case of the chiral -wave state () and . The vertical dotted line at corresponds to the grazing trajectory (for which ).

iii.2 representation

For the chiral -wave state corresponding to the IREP , we obtain from Table 1 the following expressions for the gap functions in the two domains: and . Therefore,

(19)

and .

. In this case, Eq. (9) takes the following form:

(20)

The ABS energy is discontinuous at and also at , see Appendix A. It has two lines of zeros at and another one in the basal plane, i.e., at . The point zeros at and correspond to the trajectories parallel to the DW, for which the Andreev approximation is not applicable. In the upper panels of Fig. 3, we show the ABS energy for and .

The DOS is given by Eq. (8). Following the same steps as in the previous subsection, we obtain a constant DOS:

(21)

for , see the bottom panel of Fig. 3. For , we have from Eq. (20): , and Eq. (8) can be reduced to the form

where . The last integral can be easily evaluated and we arrive at the following final expression for the DOS:

(22)

which diverges logarithmically at , as shown Fig. 3. This divergence is nothing but the van Hove singularity due to the saddle points in the ABS dispersion at and , i.e., for perpendicular to the DW.

. We obtain from Eq. (11):

(23)

The energy is discontinuous at and , see Appendix A. As for the IREP, the ABS dispersion as a function of simply shifts upon changing , therefore the DOS does not depend on . A straightforward calculation yields the following result:

(24)

The ABS energy for and the DOS for any are shown in Fig. 4.

Figure 3: (Color online) The ABS energy as a function of the direction of semiclassical propagation (upper panels) and the corresponding DOS (lower panels), for and , in the case of the chiral -wave state () and . The grazing trajectories (for which ) are shown at and by horizontal dotted lines.
Figure 4: (Color online) The ABS energy as a function of the direction of semiclassical propagation for (upper panel) and the corresponding DOS (lower panel), in the case of the chiral -wave state () and . The vertical dotted line at corresponds to the grazing trajectory (for which ).

iii.3 representation

For the chiral -wave state corresponding to the IREP , we obtain from Table 1 the following expressions for the gap functions in the two domains: and . Therefore,

(25)

and .

. In this case, Eq. (9) takes the following form:

(26)

This expression has discontinuities at and also at , see Appendix A. It has four lines of zeros at and , and another one at . The isolated second-order point zeros at and correspond to the trajectories parallel to the DW, for which the Andreev approximation is not applicable. The ABS energy for and is shown in the upper panels of Fig. 5.

The DOS for , see Eq. (8), can be reduced to the following form:

(27)

where is the Heaviside step function and . Since the function attains its maximum at , the DOS vanishes at . The integral in Eq. (27) is evaluated numerically. The logarithmic van Hove singularity in the DOS at is due to the saddle points in the ABS dispersion in the basal plane, at and .

One can easily show that the DOS has a zero-energy singularity at all values of . Indeed, we have

away from the spectrum discontinuities. This last expression has the following zeros: (i) , whose contribution to the DOS is nonsingular, due to the factor in front of the -function in Eq. (8); (ii) and , which corresponds to a maximum (minimum) of ; and (iii) and , which corresponds to the saddle points of . It is the saddle points, which are located at the four perpendicular directions in the basal plane where the lines of zeros of intersect, that produce the van Hove singularity at . The DOS for and are shown in Fig. 5.

. It follows from Eq. (11) that

(28)

which is discontinuous at and , see Appendix A. The ABS dispersion as a function of shifts upon changing , therefore, the DOS does not depend on and we obtain from Eq. (10):

(29)

where . The integral here is calculated numerically. The ABS energy for and the DOS for any are shown in Fig. 6. The logarithmic singularity in the DOS at comes from the saddle points in the ABS dispersion at , i.e., for perpendicular to the DW.

Figure 5: The ABS energy as a function of the direction of semiclassical propagation (upper panels) and the corresponding DOS (lower panels), for and , in the case of the chiral -wave state () and . The grazing trajectories (for which ) are shown at and by horizontal dotted lines.
Figure 6: (Color online) The ABS energy as a function of the direction of semiclassical propagation for (upper panel) and the corresponding DOS (lower panel), in the case of the chiral -wave state () and . The vertical dotted line at corresponds to the grazing trajectory (for which ).

iii.4 representation

Finally, we consider the chiral -wave state corresponding to the IREP , in which case and . Therefore,

(30)

and .

. In this case, Eq. (9) takes the following form:

(31)

which is discontinuous at and also at , see Appendix A. It has four lines of zeros at and . The isolated second-order point zeros at and correspond to the trajectories parallel to the DW, for which the Andreev approximation is not applicable. The ABS energy for and is shown in the upper panels of Fig. 7.

The DOS, see Eq. (8), can be calculated analytically for and . Following the same steps as in the previous subsections, we obtain:

(32)

for , and

(33)

for , see Fig. 7.

. We obtain from Eq. (11):

(34)

The discontinuities of the ABS energy are located at and also at , see Appendix A. As in the previous subsections, the ABS dispersion as a function of shifts upon changing , the DOS does not depend on , and we obtain:

(35)

where . The ABS energy for and the DOS for any are shown in Fig. 4. The logarithmic singularity in the DOS at comes from the saddle points in the ABS dispersion at , i.e., for perpendicular to the DW.

Figure 7: The ABS energy as a function of the direction of semiclassical propagation (upper panels) and the corresponding DOS (lower panels), for and , in the case of the chiral -wave state () and . The grazing trajectories (for which ) are shown at and by horizontal dotted lines.
Figure 8: (Color online) The ABS energy as a function of the direction of semiclassical propagation for (upper panel) and the corresponding DOS (lower panel), in the case of the chiral -wave state () and . The vertical dotted line at corresponds to the grazing trajectory (for which ).

iii.5 Topological origin of the ABS zero modes

The number of zero-energy ABS localized at the DW separating degenerate chiral states is determined by the difference between topological invariants characterizing the superconducting states in the bulk of the domains, which is known as the bulk-boundary correspondence.Volovik () As an illustration of this statement, we focus on the case of . To define the appropriate topological invariant, we introduce the Matsubara-like Green’s function of the Bogoliubov quasiparticles in the bulk:

(36)

where is imaginary “frequency”, takes values in the 3D Brillouin zone, and

is the Bogoliubov-de Gennes (BdG) Hamiltonian, with the gap function . Since we consider only the singlet pairing and the triplet pairing with , see Sec. II, the spin channels are decoupled and the BdG equations are reduced to a two-component (electron-hole, or Nambu) form for each spin. The BdG Hamiltonian can be written in the form , where are the Pauli matrices in the Nambu space and

The eigenvalues of are given by , where the energy of the Bogoliubov fermionic excitations.

At given , regarded as a parameter, one can define the following topological invariant:Volovik ()

(37)

Here “tr” stands for the Nambu matrix trace, the powers of the 1-form should be understood in the sense of combined exterior and matrix multiplication, and the integration is performed over and , with taking values in the 2D cross-section of the Brillouin zone by the constant plane. After some algebra, we obtain:

(38)

where . We assume that the superconducting pairing is BCS-like and effective only near the Fermi surface. At given , this results in the gap function being nonzero only near the Fermi line , which is the intersection of the Fermi surface and the constant plane. We represent the gap function in the form , and assume that there are no gap nodes, i.e. the gap magnitude does not vanish anywhere on the Fermi line. Then it follows from Eq. (38) that

(39)

therefore the topological invariant (37) is nothing but the phase winding number of the gap function around the cross-section of the Fermi surface at given . Assuming a spherical Fermi surface, the cross-section is a circle of radius . For the superconducting states considered above, the topological invariant (39) takes opposite nonzero values for the states of opposite chirality, see Table 2. The topological invariants are not defined at the bulk gap nodes, i.e. at for all four IREPs and additionally at for the IREPs and .

According to Eqs. (13), (20), (26), and (31), the ABS energy for all four IREPs can be written in the following form:

(40)

where for and for , and the function depends on the IREP. At fixed , the last expression vanishes at some values of , corresponding to the ABS zero modes. It is easy to see that there are zero modes: for they correspond to , while for they correspond to and .

One can define the algebraic number of the ABS zero modes as the number of positive-velocity modes minus the number of negative-velocity modes. According to the bulk-boundary correspondence (Ref. Volovik, ), is equal to the difference between the topological invariants in the bulk of the two domains:

(41)

Expressing the ABS energy, see Eq. (40), in terms of , one can show that the ABS zero modes propagate along the DW in the same direction: , therefore . On the other hand, it follows from Table 2 that and , which means that Eq. (41) is indeed satisfied. Taking into account the doubling of the degrees of freedom due to spin, the total number of the ABS zero modes localized near the DW is equal to , at given . Note that the same topological argument can be used to prove the existence of zero-energy ABS near the surface of an unconventional superconductor, see Ref. ZBCP, . For UPt, it was done recently in Ref. GN15, .

IREP gap function
Table 2: Topological invariant, Eq. (39), for the chiral states corresponding to the 2D IREPs of .

Iv Conclusion

We have found that the DWs separating degenerate TRS-breaking superconducting states in a 3D hexagonal crystal always create the quasiparticle ABS, for all directions of the semiclassical propagation. We have considered all four 2D IREPs of the point group (two singlet and two triplet cases) and two orientations of the DW, parallel and perpendicular to the axis. The ABS spectrum strongly depends on the order parameter symmetry and the DW orientation. Additionally, it is affected by the Josephson phase difference between the domains, which is determined by the microscopic parameters. If the DW is parallel to the axis, then there is a significant difference between the chiral states (corresponding to ) and (corresponding to ), which can be treated analytically. The spectrum of the DW ABS’s can be probed in tunneling experiments by measuring their DOS, which has very different energy dependence from that of the bulk quasiparticles. We have calculated the DOS per unit area of the DW and found a widely varying behaviour, the most prominent feature being the logarithmic van Hove singularity at zero energy, which is present in several cases.

Despite the qualitative sensitivity of the DOS to the microscopic parameters that cannot be easily controlled in experiment, we can still make some firm predictions for the DW effects on the tunneling measurements in UPt. First, there is strong evidence that the gap symmetry in the phase of UPt is described by the chiral -wave state corresponding to the IREP . If this is the case, then our results in Sec. III.3 indicate that the zero-energy singularity in the DOS is a universal feature, which, in contrast to the other three IREPs, is present for both orientations of the DW and for all values of . Second, if the DW is perpendicular to the axis, then the DOS does not actually depend on , showing different behaviour for the four IREPs: a broad dome-like maximum for , a constant for , the zero-energy singularity with two sharp edges for , and the zero-energy singularity without sharp edges for . We hope that these features can be directly probed in tunneling experiments, thus shedding light on the presence of the DWs as well as the underlying pairing symmetry in UPt and other hexagonal superconductors.

Acknowledgements.
This work was supported by a Discovery Grant from the Natural Sciences and Engineering Research Council of Canada.

Appendix A Discontinuities of the ABS spectrum

In all cases studied in this paper the gap function has the same magnitude on both sides of the DW: