Fieldangledependent lowenergy excitations around a vortex in the superconducting topological insulator CuxBi2Se3
\abstWe study the quasiparticle excitations around a single vortex in the superconducting topological insulator CuBiSe, focusing on a superconducting state with point nodes. Inspired by the recent Knight shift measurements, we propose two ways to detect the positions of point nodes, using an explicit formula of the density of states with KramerPesch approximation in the quasiclassical treatment. The zeroenergy local density of states around a vortex parallel to the axis has a twofold shape and splits along the nodal direction with increasing energy; these behaviors can be detected by the scanning tunneling microscopy. An angular dependence of the density of states with a rotating magnetic field on the  plane has deep minima when the magnetic field is parallel to the directions of point nodes, which can be detected by angularresolved heat capacity and thermal conductivity measurements. All the theoretical predictions are detectable via standard experimental techniques in magnetic fields. The discovery of topological insulators and superconductors attracts a great deal of attention in condensed matter physics. A nontrivial topological feature of the bulk state shows up in the edge boundaries with closing of insulating or superconducting gaps.[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11] The Majorana fermion (a gapless zeroenergy quasiparticle) in topological superconductors is a curious particle whose annihilation and creation operators are identical. A quest for the Majorana fermions is now an exciting issue.
A typical topological insulator BiSe characterized by a topological invariant shows superconductivity with Cu intercalating[12, 13, 14, 15]. The copperdoped bismuthselenium compounds are candidates for bulk topological superconductors. Although their properties are studied actively,[16, 17, 18] identifying the gapfunction type is an unsettled issue. Several groups observed zerobias conductance peaks (ZBCPs) by point contact spectroscopy[16, 17]. The measured ZBCPs can be evidence of topological superconductivity, since the topologically protected gapless Majorana fermions form at the edges. However, Levy et al. reported scanning tunneling spectroscopy (STS) measurements in CuBiSe. The tunneling spectrum shows that the density of states (DOS) at the Fermi level is fully gapped, without any ingap states; a fullygapped nontrivial superconductivity occurs in CuBiSe. Thus, clear evidence of the topological superconductivity is now in great demand.
The recent Knightshift measurements[19] show the presence of inplane anisotropy in CuBiSe; this anisotropy can be related to the nodal characters of the superconducting state on  plane[20], since the electronic structure in normal states is isotropic. A theoretical model[15, 21] of this superconductor predicts that within onsite pairing interaction there are four different pairing symmetries, classified by the pointgroup representations (, , , and ). The oddparity state (socalled ) is fully gapped, whereas the oddparity and states ( and , respectively) have the point nodes. Thus, the anisotropic Knight shift suggests that CuBiSe is a topological superconductor with breaking the rotational isotropy.
An anisotropy feature originating from nodes in a superconducting state is detectable in CuBiSe by various ways. The angle dependence of the thermal conductivity in the basal  plane predicts[22] distinct strong anisotropy in the representation, without magnetic fields. Detecting quasiparticle excitations around a vortex is also useful for determining point nodes. A pattern of the local density of states (LDOS) around a vortex, which can be detected by scanning tunneling microscopy/spectroscopy(STM/STS), have strong anisotropy in unconventional superconductors, since the spectrum due to the Andreev bound states around a vortex is sensitive to the anisotropy of the pair potential[26, 27, 28, 23, 24, 25]. Furthermore, thermal transport measurements such as heat capacity and thermal conductivity measurements with a rotating magnetic field are probes for node positions in unconventional superconductors[29, 30].
Let us here discuss what is an efficient and intuitive way to seek the unconventional properties of the topological superconductor. A clue is the quasiclassical approach of superconductivity[31]. The author showed that the oddprity superconductivity in topological superconductors turns to the spintriplet one in terms of the quasiclassical treatment. The original massive Dirac Bogoliubovde Gennes (BdG) Hamiltonian derived from a tightbinding model represented by an matrix is reduced to a matrix. This indicates that lowenergy nontrivial quasiparticle excitations in the topological superconductors can be analyzed using established theoretical techniques of the spintriplet superconductors. With the use of the Riccati parametrization, equations of motion of quasiclassical Green’s function represented by a matrix becomes two Riccatitype differential equations represented by matrices. In addition, if the effective gap function is unitary (), the quasiclassical Green’s function is obtained by solving two Riccatitype differential equations with scalar coefficients. We use KramerPesch approximation (KPA) for these scalar Riccati equations to analyze the fieldangledependent lowenergy excitations around a vortex, which has been used in various kinds of unconventional superconductors[25, 33, 32].
In this paper, we study the quasicparticle excitations around a vortex in the threedimensional superconducting topological insulator CuBiSe, focusing on a superconducting state with point nodes. Using the quasiclassical treatment, we derive Riccatitype firstorder differential equations with scalar coefficients related to the original massive Dirac Hamiltonian. The KPA to study the Andreev bound states leads to an explicit formula of the LDOS around a single vortex. We propose two ways to detect the inplane anisotropy of the topological superconductivity. In the magnetic field parallel to the axis, a twofold shape of the LDOS pattern around a vortex at zero energy splits along the nodal direction with increasing energy; this is detectable in the standard STM measurements. In the magnetic field perpendicular to the axis, an angular dependence of the DOS on a rotating magnetic field on the  plane has deep minima when the magnetic field is parallel to the directions of point nodes; the angularresolved heatcapacity and thermalconductivity measurements may show this behavior.
We study a model of a topological superconductor in a threedimensional system, with meanfield approximation. As for the normal parts, the effective Hamiltonian with a strong spinorbit coupling around the point is equivalent to that of the massive Dirac Hamiltonian with the negative Wilson mass term. The Bogoliubovde Gennes (BdG) meanfield Hamiltonian of the superconductivity is , with the 8component column vector . The corresponding BdG equations with matrix eigenequations are[22, 31, 34, 35, 36]
(1) 
with the normalstate Hamiltonian
(2) 
and the pairing potential matrix . The matrices () are the Gamma matrices in the Dirac basis. Using the Pauli matrices and in orbital and spin spaces, we have and , with a unit matrix . Within onsite interaction, the pair potential must fulfill the relation owing to the fermionic property. We have six possible gap functions classified by a Lorentztransformation property[22]; they are classified into a pseudoscalar, a scalar, and a polar vector (fourvector). A direction of point nodes is parallel to that of polar vectors, since the rotational symmetry is preserved around this direction (See, Appendix B in Ref. \citenNagaiThermal). Thus, we propose a pair potential with point nodes in direction on  plane,
(3)  
(4) 
Diagonalizing the BdG Hamiltonian in the uniform system, we find that the pointnodes are located at
(5) 
The momentumindependent gap function makes anisotropic energy spectrum, although the Fermi surface of the normal state is spherical in this model. In this paper, we focus on physics originating from the point nodes in the topological superconductor, from a phenomenological point of view, although the microscopic superconducting mechanism supporting the Knightshift measurements is quite interesting.
Here, we use the quasiclassical approach, to reduce the matrix to a one[31]. Using this approach, we show the fact that an oddparity superconductivity in topological superconductors maps into a spintriplet one. According to Ref. \citenNagaiQuasi, the quasiclassical Andreev equations with matrix eigenequations are
(6) 
with the effective pairing potential matrix . Considering the oddparity gap functions, is defined by . The odd parity gap function with point nodes shown in Eq. (4) is represented by , with the Fermi velocity . The corresponding quasiclassical Green’s function is
(7)  
(8) 
where the matrix coefficients and obey matrix Riccati equations, [25, 37, 38, 39, 40, 41, 42]
(9)  
(10) 
Here, denotes the Fermion matsubara frequency. Considering clean superconductors in the typeII limit, we neglect the selfenergy part of the Green’s function and the vector potential. The following results do not change qualitatively, with taking impurity selfenergies, as far as the pointnodes do not vanish in a dirty system. Assuming that the effective pairing potential matrix is a unitary matrix (), we obtain scalar Riccati equations,
(11)  
(12) 
where the scalars , and are defined by , and , respectively. Since these equations contain only through , the above equations reduces to a onedimensional problem on a straight line, the direction of which is given by that of the Fermi velocity . Considering the cylindrical coordinate frame with in the real space, the pair potential does not depend on in the Riccati equations, and hence the Riccati equations can be rewritten as
(13)  
(14) 
where is an amplitude of the vector perpendicular to the axis taken by projecting the Fermi velocity. Here, the argument is along the direction parallel (perpendicular) to .
Let us introduce the KPA as an efficient method to analyze the angleresolved experiments. The zeroenergy density of states around a vortex given by KPA is consistent with that of direct numerical calculations, as seen, e.g., in Ref. \citenNagaiPRL. We expand the Riccati equations up to first order with respect to energy and the imaginary part of the pair function[43, 44]. For a single vortex, we adopt the spatial variation of the pair potential expressed by
(15) 
with the cylindrical coordinate . The LDOS with the use of KPA is given as
(16) 
where is the smearing factor, is the Fermisurface area element, and
(17)  
(18)  
(19) 
The quantity is defined by
(20) 
where the function is the modified Bessel function of the second kind.
Now, we study the LDOS pattern around a vortex in the magnetic field parallel to the axis. For simplicity, we adopt Polarvector type gap function (socalled in Ref. \citenPhysRevLett.107.217001), where the point nodes are located on the axis (). We set the smearing factor Near zero energy, the bound states spread to the antinodal direction as shown in Fig. 1(a) and Fig. 1(b). In contrast, with increasing energy, the bound states spread to the nodal direction as shown in Fig. 1(c) and Fig. 1(d). These results indicate that a twofold shape of the LDOS pattern at zero energy splits along the nodal direction with increasing energy. In Figs.1(a)(b), the shape is elliptic spreading in the direction. On the other hand, in Figs.1(c)(d), the peak position is located on the axis. The shape of the LDOS pattern was discussed in terms of the enveloping curve of the quasiparticle paths[24]. The LDOS in wave superconductors spreads to the nodal direction at the zero energy, since there is the asymptote of the trajectory in the nodal direction[24]. However, the trajectory in this system is a parabola. The zeroenergy LDOS spreads to the antinodal direction, since the antinodal quasiparticle paths are dense near the vortex because of the sharp curve of trajectory. The split of the LDOS pattern is evidence of unconventional superconductivity, since they have been observed in the superconductors with the anisotropic superconducting gap structure such as NbSe[26, 45] and YNiBC[27, 33].

Next, we examine the angularresolved DOS in the magnetic field perpendicular to the axis. The DOS is calculated by
(21) 
Here, is the realspace average around a vortex where [ and ]. We set the spatial cutoff length , which is comparable to the neighboring vortex distance as the magnetic field .

The zeroenergy DOS drastically decreases, when the magnetic field is applied in the parallel direction to the direction of the nodal quasiparticles (i.e., ), as shown in Fig. 2. Here, denotes the direction of magnetic fields. The order of this reduction is similar to the case of wave superconductivity with a cylindrical Fermi surface[32]. In terms of the Doppler shift method, which has been frequently used to analyze experiments[29, 30], the quasiparticles around a vortex have the Dopplershift energy , where means the superfluid velocity around a vortex and . The quasiparticles near the gap nodes mainly contribute to the DOS, since they are excited only by the Dopplershift energy. In the magnetic field parallel to the gapnode direction, the nodal quasiparticles have no Dopplershift energy, because . In terms of the KPA, the DOS originating from the bound states around a vortex reduces in the magnetic field parallel to the direction of the nodal quasiparticles, since the quasiparticles along the magnetic field are not bound. The angularresolved DOS in Fig. 2 can be observed by angularresolved heat capacity and thermal conductivity measurements[29, 30], in topological superconductors with pointnodes. In this paper, we do not take the contributions in the DOS originating from the surface bound states into account. Their contributions may be predominant in a small single crystal; indeed, the surface bound states lead to strong anisotropy in the thermal conductivity without a magnetic field, as seen in Ref. \citenNagaiThermal. Our results indicate that the anisotropy from the point nodes can be detected by the angularresolved DOS, even in the bulk, with a rotating magnetic field on the  plane.
In conclusion, we proposed two ways to detect the position of pointnodes in the model of superconducting topological insulator CuBiSe with the quasiclassical approach of superconductivity. We reduced matrix eigenequations to two scalar Riccatitype differential equations. Furthermore, we applied the KPA, to obtain an explicit formula of the DOS. A twofold shape of the LDOS pattern around a vortex at zero energy splits along the nodal direction with increasing energy in the magnetic fields parallel to the axis. The DOS with a rotating magnetic field on the  plane has deep minima in the magnetic field parallel to the directions of point nodes. All the theoretical predictions are detectable via standard experimental techniques, including the STM and the anglerresolved thermal transport measurements.
We thank Y. Ota, H. Nakamura and M. Machida for helpful discussions and comments. The calculations were performed using the supercomputing system PRIMERGY BX900 at the Japan Atomic Energy Agency. This study was supported by JSPS KAKENHI Grant Number 26800197.
References
 B. A. Bernevig, T. L. Hughes, and S.C. Zhang: Science 314 (2006) 1757.
 Y. L. Chen, Z. K. Liu, J. G. Analytis, J.H. Chu, H. J. Zhang, B. H. Yan, S.K. Mo, R. G. Moore, D. H. Lu, I. R. Fisher, S. C. Zhang, Z. Hussain, and Z.X. Shen: Phys. Rev. Lett. 105 (2010) 266401.
 L. Fu and C. L. Kane: Phys. Rev. B 76 (2007) 045302.
 L. Fu, C. L. Kane, and E. J. Mele: Phys. Rev. Lett. 98 (2007) 106803.
 M. Z. Hasan and C. L. Kane: Rev. Mod. Phys. 82 (2010) 3045, and references therein.
 C. L. Kane and E. J. Mele: Phys. Rev. Lett. 95 (2005) 146802.
 M. König, S. Wiedmann, C. Brüne, A. Roth, H. Buhmann, L. W. Molenkamp, X.L. Qi, and S.C. Zhang: Science 318 (2007) 766.
 K. Kuroda, M. Ye, A. Kimura, S. V. Eremeev, E. E. Krasovskii, E. V. Chulkov, Y. Ueda, K. Miyamoto, T. Okuda, K. Shimada, H. Namatame, and M. Taniguchi: Phys. Rev. Lett. 105 (2010) 146801.
 J. E. Moore and L. Balents: Phys. Rev. B 75 (2007) 121306.
 A. Nishide, A. A. Taskin, Y. Takeichi, T. Okuda, A. Kakizaki, T. Hirahara, K. Nakatsuji, F. Komori, Y. Ando, and I. Matsuda: Phys. Rev. B 81 (2010) 041309.
 T. Sato, K. Segawa, H. Guo, K. Sugawara, S. Souma, T. Takahashi, and Y. Ando: Phys. Rev. Lett. 105 (2010) 136802.
 Y. S. Hor, A. J. Williams, J. G. Checkelsky, P. Roushan, J. Seo, Q. Xu, H. W. Zandbergen, A. Yazdani, N. P. Ong, and R. J. Cava, Phys. Rev. Lett. 104, 057001 (2010).
 L. A. Wray, S.Y. Xu, Y. Xia, Y. S. Hor, D. Qian, A. V. Fedorov, H. Lin, A. Bansil, R. J. Cava, and M. Z. Hasan, Nat Phys 6, 855 (2010).
 M. Kriener, K. Segawa, Z. Ren, S. Sasaki, and Y. Ando, Phys. Rev. Lett. 106, 127004 (2011).
 L. Fu and E. Berg, Phys. Rev. Lett. 105, 097001 (2010).
 S. Sasaki, M. Kriener, K. Segawa, K. Yada, Y. Tanaka, M. Sato, and Y. Ando, Phys. Rev. Lett. 107, 217001 (2011).
 T. Kirzhner, E. Lahoud, K. B. Chaska, Z. Salman, and A. Kanigel, Phys. Rev. B 86, 064517 (2012).
 N. Levy, T. Zhang, J. Ha, F. Sharifi, A. A. Talin, Y. Kuk, and J. A. Stroscio, Phys. Rev. Lett. 110, 117001 (2013).
 K. Matano, H. Yamamoto, K. Ueshima, F. Iwase, M. Kriener, K. Segawa, Y. Ando, G.q. Zheng, JPS 2014 Annual (69th) Meeting, 28aBF12.
 T. Hashimoto, K. Yada, A. Yamakage, M. Sato and Y. Tanaka, J. Phys. Soc. Jpn. bf 82, 044704 (2012).
 Y. Nagai, H. Nakamura, and M. Machida, arXiv:1211.0125v2 to be published in J. Phys. Soc. Jpn.
 Y. Nagai, H. Nakamura, and M. Machida, Phys. Rev. B 86, 094507 (2012).
 N. Hayashi, M. Ichioka and K. Machida, Phys. Rev. B 56, 9052 (1997).
 Y. Nagai, Y. Ueno, Y. Kato and N. Hayashi, J. Phys. Soc. Jpn. 75, 104701 (2006).
 Y. Nagai, Y. Kato and N. Hayashi, J. Phys. Soc. Jpn. 75, 043706 (2006).
 H. F. Hess, R. B. Robinson and J. V. Waszczak, Phys. Rev. Lett. 64, 2711 (1990).
 H. Nishimori, K. Uchiyama, S. Kaneko, A. Tokura, H. Takeya, K. Hirata and N. Nishida, J. Phys. Soc. Jpn. 73 3247 (2004).
 F. Gygi and M. Schluter, Phys. Rev. B 43, 7609 (1991).
 T. Sakakibara, A. Yamada, J. Custers, K. Yano, T. Tayama, H. Aoki and K. Machida, J. Phys. Soc. Jpn. 76, 051004 (2007).
 Y. Matsuda, K. Izawa, and I. Vekhter, J. Phys. Condens. Matter 18, R705 (2006).
 Y. Nagai, H. Nakamura, and M. Machida, J. Phys. Soc. Jpn. 83, 053705 (2014).
 Y. Nagai, and N. Hayashi, Phys. Rev. Lett. 101, 097001 (2008).
 Y. Nagai, Y. Kato, N. Hayashi, K. Yamauchi, and H. Harima, Phys. Rev. B 76, 214514 (2007).
 L. Hao and T. K. Lee, Phys. Rev. B 83, 134516 (2011).
 A. Yamakage, K. Yada, M. Sato, and Y. Tanaka, Phys. Rev. B 85, 180509 (2012).
 T. Mizushima, A. Yamakage, M. Sato, and Y. Tanaka, arXiv:1311.2768 (unpublished).
 Y. Kato, J. Phys. Soc. Jpn. 69, 3378 (2000).
 Y. Nagato, K. Nagai, and J. hara, J. Low Temp. Phys. 93, 33 (1993).
 S. Higashitani and K. Nagai, J. Phys. Soc. Jpn. 64, 549 (1995).
 Y. Nagato, S. Higashitani, K. Yamada, and K. Nagai, J. Low Temp. Phys. 103, 1 (1996).
 N. Schopohl and K. Maki, Phys. Rev. B 52, 490 (1995).
 N. Schopohl, arXiv:condmat/9804064 (unpublished).
 A. S. Melnikov, D. A. Ryzhov, and M. A. Silaev, Phys. Rev. B 78, 064513 (2008).
 Y. Nagai, H. Nakamura, and M. Machida, Phys. Rev. B 83, 104523 (2011).
 N Hayashi, M Ichioka, and K Machida, Phys. Rev. Lett. 77, 4074 (1996).