Vortices in chiral, spin-triplet superconductors and superfluids

Superconductors exhibit unconventional electronic and magnetic properties if the Cooper pair wave function breaks additional symmetries of the normal phase. Rotational symmetries in spin- and orbital spaces, as well as discrete symmetries such as space and time inversion, may be spontaneously broken. When this occurs in conjunction with broken global U(1) gauge symmetry, new physical phenomena are exhibited below the superconducting transition that are characteristic of the broken symmetries of the pair condensate. This is particularly true of vortices and related defects. Superconductors with a multi-component order parameter exhibit a variety of different vortex structures and closely related defects that are not possible in condensates belonging to a one-dimensional representation. In this article we discuss the structure of vortices in Fermionic superfluids and superconductors which break chiral symmetry, i.e. combined broken time-inversion and 2D parity. In particular, we consider the structure of vortices and defects that might be realized in thin films of 3He-A and the layered superconductor Sr2RuO4, and identify some of the characteristic signatures of broken chiral symmetry that should be revealed by these defects.


Introduction
The BCS theory of superconductivity [1], combined with refinements over several decades [2,3,4,5], ranks among the major achievements in theoretical physics during the last century. The central feature of BCS theory is Cooper pair condensation -i.e. the macroscopic occupation of a single quantum state of bound pairs of fermions. In its simplest form for s-wave, spin-singlet pairing the "order parameter", or Cooper pair amplitude, is given a complex scalar function of the relative coordinate, r, and center of mass position, R, of the pairs, The Cooper pair amplitude is both a measure of the condensate density, |Ψ| 2 ∼ O(N/V ), and a reflection of the spontaneously broken symmetry of the ordered phase, which in this case is the global U(1) gauge symmetry generated by Fermion particle number [6]. For homogeneous states the phase, ϑ, is the signature of the broken U(1) symmetry, i.e. a macroscopic fraction of Fermions condense into a two-particle state with the same quantum phase. The broken global U(1) symmetry also implies a family of degenerate ground states related to one another by the elements of the symmetry group. Long-wavelength spatial variations of the phase form a branch of low-lying collective excitations, the Anderson-Bogoliubov phase mode [7,8], with an acoustic dispersion relation. This is a classic example of a Goldstone mode that accompanies a broken continuous symmetry. In superconductors, Cooper pairs are charged particles with charge 2e, and the requirement of local U(1) gauge symmetry leads to a massive Goldstone boson via the Higgs mechanism, that corresponds to the longitudinal oscillation of the electromagnetic vector potential, and is a high energy plasmon excitation ‡. In neutral superfluids the phase mode is observable as fourth sound, a collective excitation of the condensate which has a velocity that is determined by the "phase stiffness", ρ s ∝ |Ψeq| 2 , or superfluid density, where Ψeq is the equilibrium condensate amplitude. This rigidity stabilizes the condensate against phase fluctuations, and is responsible for the Josephson effect [9], persistent currents [10], and flux quantization in superconductors. Indeed the standard signature of BCS condensation is the quantization of flux in units of Φ 0 = hc/2e, or the quantization of circulation in units of κ 0 = h/2m in a neutral pair condensate of Fermions of mass m. Quantized flux and quantized circulation are consequences of a fundamental constraint on the phase of the pair amplitude, where N C = 0 ± 1, ±2, . . . is the winding number of the phase around a closed contour C within the condensate. This quantization condition reflects the physical requirement that the order parameter be a single-valued field, which for a complex scalar pair amplitude is simply that the phase return to its value modulo 2π. For a multiply connected geometry this constraint and the phase stiffness implies energy barriers separating states of different winding number, and thus to quantized persistent currents. In a simply connected geometry states with N C = 0 force the condensate to effectively become multiply connected as there is necessarily one or more phase singularities interior to the contour C. These singularities, points in two dimensions or lines in three dimensions, form a spectrum of topologically stable "defects" of the order parameter, i.e. phase vortices, labelled by their winding number. These topologically stable defects are often energetically stable, or metastable because of energy barriers separating states with different winding numbers. The local structure of these defects, the vortex cores, are determined in part by their topology -their winding number in this simple case, the local spectrum of fermionic excitations and their interactions, as well as scattering of the Fermions by impurities or other sources of disorder [11]. BCS superconductors exhibit unconventional properties if the order parameter breaks additional symmetries of the normal phase. Spin-and orbital rotational invariance, as well as space and time inversion, may be spontaneously broken. When this occurs in conjunction with broken U(1) gauge symmetry, new collective excitations of the pair condensate [12,13], novel heat and transport properties [14,15], unconventional vortices [16,17] as well as point defects [18,19] characteristic of the complex symmetry breaking are possible [20]. In this article we discuss the structure of vortices in superfluid 3 He and superconductors which break spin and orbital rotation symmetries, parity and/or time-inversion symmetry. This class of pairing states is believed to describe the ground states of thin films of 3 He-A [21], the low-temperature superconducting phase of UPt 3 [22], and the layered superconductor Sr 2 RuO 4 [23].

Theoretical Formalism
The superconducting order parameter can be defined in terms of the condensate amplitude for Cooper pairs, where ψ α (r) is the fermion field operator for spin projection α. This amplitude is directly related to the anomalous Matsubara propagator [24], expressed here in terms of the center-of-mass, R, and relative (orbital) coordinate, r, of the pair, and the (imaginary) time difference, τ = τ 1 − τ 2 . T τ is the (imaginary) time-ordering operation for Fermions. Fourier transforming with respect to the orbital coordinate and time gives the anomalous propagator for pairs with orbital momentum p, center of mass position R, and Matsubara energy, ǫ n = (2n + 1)π/β for n = 0, ±1, . . ., where β = 1/k B T . The local spectral function for the pairs is obtained by analytic continuation to the real axis, iǫ n → ǫ + i0 + [24]. The orbital radius of a Cooper pair is of order ξ 0 = v f /2πk B T c , which is typically large compared with atomic scales, e.g. the Fermi wavelength, ξ 0 ≫ /p f . In this limit we can factor out the short wavelength variations of the pair amplitude and compute directly the spatial variations of the pair amplitude and order parameter on mesoscopic length scales defined by the coherence length, ξ 0 , the mean-free path, ℓ and/or the London penetration depth, λ. This separation occurs because the high-energy, ǫ ∼ E f , short-wavelength, |p| ∼ p f , properties of the Fermi liquid are unaffected by the longwavelength, low-energy pairing correlations to leading order in /p f ξ 0 , /p f ℓ, k B T c /E f etc. Thus, we define the amplitude for Cooper pairs near the Fermi surface by integrating over the low-energy band of states near the Fermi level [26], where p f is the Fermi momentum, ξ p = v f (|p| − p f ) and Ω c ≪ E f is the quasiparticle bandwidth. A detailed discussion of this procedure is given in Refs. [25,26,27,28,29]. A synopsis of the quasiclassical transport theory for equilibrium states of inhomogeneous superconductors, as well as the methods used to calculate the order parameter and electronic structure of superconducting states associated with impurities, interfaces and vortices is described below. The propagator for Cooper pairs with relative momenta on the Fermi surface given in Eq. 6 is coupled to the low-energy propagator for fermionic quasiparticles of the superconductor. The coupled equations for quasiparticles and Cooper pairs are formulated in terms of a 4× 4 matrix propagator in the combined particle-hole (Nambu) and spin space, which incorporates the the spin correlations that are required by the pairing correlations, is the four-component (Nambu spinor) Fermion field operator, r 1 = R + r/2, r 2 = R − r/2. Note also that for imaginary time evolution, ψ(r, τ ) ≡ ψ † (r, −τ ). The quasiclassical matrix propagator is then defined by integrating the over the low-energy band of states near the Fermi surface, where τ 3 is the third Pauli matrix in Nambu space, and we have renormalized by dividing by the spectral weight, a, of the normal-state quasiparticle resonance, 0 < a < 1. Its structure in Nambu space is represented by the quasiparticle propagators,ĝ andĝ, and the pair propagators,f andf. Each of these propagators is a 2 × 2 matrix in spin space denoted by small hats. The pair propagator naturally separates into spin-singlet (anti-symmetric) and spin-triplet (symmetric) amplitudes, while the natural description for the quasiparticle propagator is in terms of spin-scalar and vector components, where σ = (σ x , σ y , σ z ) are the Pauli spin matrices. The description in terms for four matrix propagators is convenient, but redundant. Fundamental symmetry relations associated with Fermion anti-symmetry and conjugation (i.e particle ↔ hole) connect the two diagonal and off-diagonal propagators, whereb tr (b † ) is the transpose (adjoint) ofb.

Transport Equation
Extensions of Landau's theory of normal Fermi liquids to include BCS pairing correlations [30], electron-phonon interactions [31], and disorder [32] culminated in the equations of Eilenberger [5] and Larkin and Ovchinnikov [33], formulated in terms of the Nambu matrix propagator, g(p f , R; ǫ n ), obeying transport-type differential equations along classical trajectories defined by the Fermi momentum, § where [. . .] is the commutator in Nambu space between the matrix propagator, g, the (imaginary) time-development operator, τ 3 ∂ τ ↔ iǫ n τ 3 , and various internal and external fields that couple to quasiparticles and pairs. Eilenberger's transport equation is obtained by an energy-and momentum-space renormalization of Dyson's equation and the full many-body equations for the self-energy functional. The normalization information contained in Dyson's equation is absent in Eilenberger's equation, but is recovered in the form of a contraint [5,33] that must be satisfied by the physically allowed solutions of the transport equation. The normalization constraint on the Matsubara propagator is This constraint can be obtained as a boundary condition which supplements the transport equation, and reflects the fact that in the absence of any disturbance the propagator reduces to that for local homogeneous equilbrium. For more detailed discussions of the normalization condition see Refs. [5,33,36]. The renormalization procedure leads to an expansion of the full many-body Green's functions and self-energies in terms of ratios of low to high energy scales, e.g. k B T c /E f , u ext /E f , or short to long wavelength scales, e.g.
/p f ξ 0 , q/p f . We use a single parameter, small, to classify the order of magnitude of various terms in the transport equation, the magnitude of the self-energy and external fields, etc. The results reported here are based on the leading order expansion of the self energy in small, particularly the mean-field pairing self-energy, ∆, that describes the inhomogeneous equilibrium state of the superconductor, and the effects disorder in terms of scattering by a dilute random distribution of impurities represented by the impurity self energy, Σ. However, we omit the Landau molecular field self-energy, which is also of order small. We do not expect the Landau molecular field to lead us to qualitatively different conclusions regarding the basic structure of vortices in chiral p-wave states, however these terms may play a role in the relative stability of various multi-vortex configurations. We also treat the mean-field pairing self energy in the weak-coupling limit. Indeed a useful operational § Extensions of the transport theory of Eilenberger to non-equilibrium superconductivity were formulated by Eliashberg [34] and by Larkin and Ovchinnikov [35]. The nonequilibrium transport theory is discussed our companion paper on vortex dynamics in this volume [40]. definition for the order parameter is the mean field pairing self energy, ∆ = 0∆ ∆ 0 (15) where∆ is the 2 × 2 spin matrix order parameter with components ∆ αβ (p f , R; ǫ n ), loosely referred to as the "gap function", is is both a measure of the pairing correlation energy and the broken symmetry of the superconducting state. The gap function is determined by the BCS self-consistency equation relating∆(p f , R; ǫ n ) to the Cooper pair amplitude,f(p f , R; ǫ n ), and the interaction,λ(p f , ǫ n , p ′ f , ǫ ′ n ), that provides the "pairing glue". In the weak-coupling limit the frequency dependence of the pairing interaction, and the pairing self-energy, is constant within the bandwidth |ǫ n | ≤ Ω c . The order parameter then satisfies the weak-coupling gap equation, where f = is the Cooper pair propagator and is the bosonic propagator responsible for the pairing glue. N f is the single-spin normal-state density of states at the Fermi level and the integration is an average over the Fermi surface normalized to d 2 p f ≡ dS p f n(p f ) = 1, where n(p f ) is the normalized angle-resolved density of states on the Fermi surface.
The impurity self-energy is a functional of the quasiclassical propagator and defined in terms of the self-consistent impurity scattering t-matrix, where n s is the density of impurity density and t(p f , p f , ǫ n ; R) is the forward scatting limit of the self-consistent t-matrix describing multiple-scattering of quasiparticles and pairs by an impurity in the superconductor, where g(p f , ǫ n ) = is the full matrix propagator and u(p f , p ′ f ), represented by the cross and dotted line, is the matrix element of the impurity potential between normal state quasiparticles with momenta p f and p ′ f on the Fermi surface. We model the disorder in terms of non-magnetic impurities that scatter quasiparticles isotropically, i.e. we retain only the s-wave contribution to the impurity potential, thus, u = u 0 τ 3 . The corresponding scattering phase shift in the s-wave channel is defined by δ 0 = tan −1 (πN f u 0 ). The t-matrix equation, and thus the impurity self-energy, can be expressed as the solution of the matrix equation, where g(p ′′ f ; ǫ n ) is the Fermi-surface average of the propagator. In the normal state, the low-energy quasiclassical propagator reduces to g N = −iπ sgn(ǫ n ) τ 3 , and thus the t-matrix descbribing the scattering of normal state quasiparticles and quasiholes near the Fermi surface becomes, In this model the disorder is characterized by the strength of impurity potential, u 0 , the mean density of scatterers, n s and the density of states for quasiparticles at the Fermi surface. The mean lifetime of a quasiparticle on the Fermi surface in the normal state, τ N , is obtained by analytic continuation of the self-energy to the real axis to obtain the retarded self-energy, The mean free path, ℓ N = v f τ N , reduces to the classic result, ℓ N = 1/n s σ, where σ is the total cross section for scattering of a quasiparticle off an impurity, σ = (4π 2 /p 2 f ) sin δ 2 0 . In terms of the dimensionless cross sectionσ ≡ sin 2 δ 0 , the weak scattering limit ("Born limit") corresponds toσ → 0, while the strong scattering ("unitarity limit") corresponds toσ = 1.
For any inhomogeneous superconducting state, the relative importance of impurity scattering is determined by two parameters, (i) a pair-breaking parameter defined by ratio of the scale of the pair correlation length and the mean free path, x = ξ 0 /ℓ N and (ii) the strength of the impurity scattering potential as measured by the dimensionless cross section,σ. For unconventional pairing states, such as chiral p-wave superconductors, disorder destroys the superconducting state at a critical level of disorder given by x c ≃ 0.28. The results reported below are based on self-consistent calculations of the order parameter, current distribution and local density of states for the clean limit and for modest levels of disorder with ℓ N = 10ξ 0 in the Born limit.

Observables
To complete the theoretical framework for our calculations of the vortex structure for Fermi-liquid superconductors and superfluids we include expressions for the equilibrium current, and local density of states in terms of the propagator. In particular, we define the normalized angle-resolved local density of states (LDOS) by dividing out the the angle-resolved normal-state DOS at the Fermi level, N f n(p f ). Thus, the normalized local spectral density due to particle-hole coherence of the pair condensate and the formation of sub-gap excitations associated with scattering and spatial inhomogeneities of the order parameter is obtained from the retarded propagators for particle and hole excitations by analytic continuation to the real axis, where the trace is over the 4×4 Nambu space. All equilibrium properties are determined by this spectrum. In particular, the equilibrium charge current associated with a topological defect in a superconductor is given by, where f (ǫ) = 1/(e βǫ + 1) is the Fermi distribution. The orbital magnetization density is then given by ∇ × m = − 1 c j(R). For equilibrium states it is convenient to calculate the total current by transforming to the Matsubara representation, as shown in Eq. 24. However, if we are interested in the contributions to the current from specific regions of the excitation spectrum then the formulation in terms of the LDOS is more useful. In particular, the angle-resolved spectral current density is defined as the net current density carried by states at ±p f and energy ǫ [11],

Ricatti Equations
The order parameter and self-energy must be determined self-consistently with the solution of the transport equation for the propagator. An efficient method for solving the transport equation is based on a parameterization for the propagator that automatically satisfies the normalization constraint, and thus eliminates the possibility of spurious solutions. The parameterization is defined by [37,38,39,21] where the prefactor is given by The coherence amplitudesγ andγ are 2 × 2 matrices in spin space which obey Ricattitype equations, and are simply related to the particle-and hole-like projections of the off-diagonal where the projection operators in Nambu space for particle-like (+) and hole-like (-) excitations are given by [36], In particular, P 2 ± = P ± and P + P − = P − P + = 0 follow directly from the normalization condition. The physical interpretation of these projectors, and the Ricatti amplitudes, is easily established by comparing the action of the projection operators on a general Nambu spinor with the particle-like and hole-like solutions of the Bogoliubov or Andreev equations. The fundamental symmetry relation in Eqs. 11-12 for the particle and hole components of the quasiclassical propagator also imply symmetry relations relating the two Ricatti amplitudes,γ The utility of the Ricatti equations is that they provide an efficient approach to solving the transport equation. The Ricatti equations are easily integrated numerically because they have numerically stable solutions. The extension of this formalism applicable to vortex dynamics is described in a companion article for this special collection, Ref. [40], as well as Refs. [38,39,41].

Calculational Methods
The calculations for topological defects reported here are obtained from self-consistent solutions for the propagator, order parameter and self energy. For two-dimensional layered superconductors or superfluids we define all relevant functions, e.g. g(p f , R, ǫ n ), on a discrete set of points in both position and momentum space. The center of mass coordiate is discretized as R = h(ma + mb), where h is the lattice spacing, a and b are lattice unit vectors and (m, n) is a pair of integers defining a discrete lattice point. A typical choice for the lattice spacing is h = 0.25 ξ 0 , and the dimensions of the largest grids used to carry out these calculations was 120 × 120 lattice sites. The larger grids were used to investigate the stability of vortices multiple units of phase winding. We also used both rectangular and hexagonal lattices. For calculations of isolated defects we place the defect initially at the origin and enforce the constraint on the global phase winding by fixing the order parameter by its asymptotic form for |R| → ∞ along an exterior contour as shown in Fig. 1. At each grid point we use a fourth order Runge-Kutta algorithm to integrate the transport equation along a discrete set of trajectories representative of the Fermi surface and the basis functions defining the pairing symmetry. For the calculations reported here we used sixteen trajectories equally spaced in momentum space around the Fermi surface, which for simplicity we assumed to be cylindrical. However it is straight forward to implement more detailed Fermi surface geometries if needed (see Ref. [42]). We also calculate the propagator for a finite set of Matsubara frequencies. The maximum Matsubara frequency required to achieve high precision increases with decreasing temperature. For T /T c = 0.4 eight Matsubara frequencies is typically sufficient. A calculation is initialized by a "seed" for the order parameter field, with the topological constraint implemented on the boundary as described above. The propagator is then calculated for all the grid points, and for the set of Fermi surface trajectories and Matusbara frequencies at each grid point. The order parameter and impurity self energy are then calculated from discretized forms of Eqs. 16, 17 and 19. The updated values of the order parameter and impurity self-energy are then used as new inputs to the transport equation, which is solved again for each grid point to obtain an improved solution for the propagator. This iterative procedure is continued until it converges to a specified precision for the deviations of the order parameter and self energy between iterations. The final results for the equilibrium order parameter and excitation spectrum are found to be insensitive to the initial order parameter field, except of course for the constraint on the phase winding.

Pairing Symmetry
Fermi statistics requires that pair amplitude, f αβ , the pairing interaction, λ αβ;γρ , and thus the gap function, ∆ αβ , obey the anti-symmetry condition, For materials in which the normal metallic state possesses inversion symmetry the pairing interaction separates into even-(g) and odd-parity (u) channels, which are respectively anti-symmetric and symmetric under spin exchange and correspond to the spin-singlet (S = 0) and spin-triplet (S = 1) pairing channels, Furthermore, the pairing interaction separates into a sum over invariant bilinear products of basis functions for the irreducible representations, Γ, of the point group, for both even-(η Γν (p f )) and odd-parity ( η Γν (p f )) sectors, These basis functions are the eigenfunctions of the linearized gap equation, while λ Γ is the eigenvalue that determines the instability temperature for Cooper pairs belonging to the the irreducible representation, Γ, i.e. 1/λ Γ = ln(1.13Ω c /T Γ c ). Except in rare cases of accidental degeneracy, or a weakly broken symmetry [46], the pairing interactions are well separated and the most attractive channel determines both the transition temperature, T c and the basis functions, {η Γν }, for the irreducible representation, Γ, that determines the pairing state(s) of the superconductor. The index ν labels the basis functions for the representation Γ with dimension d Γ . Although the exact basis functions depend on details of the anisotropic Fermi surface and pairing interaction, representative basis functions which exhibit all the broken symmetry properties of the irreducible representation are easily constructed [47]. Below we list the irreducible representations and representative basis functions for the group D 4h , appropriate for Sr 2 RuO 4 in Table 2.4.
Unconventional pairing states which have no equal-time average, i.e. 'odd-frequency' pairing states, obey a more general anti-symmetry condition which includes the imaginary time coordinate [43]. Such states have the parity assignments for total spin quantum numbers interchanged. See Berezinskii [44], and for a more recent discussion of these states, Balatsky and Abrahams [45]. Thus, barring accidental near degeneracy [48], the order parameter, ∆ αβ (p), is defined by a single representation, and is either even-or odd-parity, and therefore spinsinglet or spin-triplet, respectively,P with ∆(p f ) = ∆(−p f ) and ∆(p f ) = − ∆(−p f ). The expansion in the eigenfunctions belonging to the dominant representation Γ (either even or odd-parity) of the linearized gap equation gives,

Odd-Parity, Spin-Triplet Pairing
The prototype for odd-parity, spin-triplet pairing is realized in the neutral Fermi superfluid 3 He. In particular the A-phase of superfluid 3 He is identified as an equal-spin pairing state with a p-wave order parameter of the form [52], where ℓ = m × n is the quantization axis for the orbital angular momentum of the pairs, i.e. ℓ · L orb (m + in) ·p = + (m + in) ·p,p = p f /|p f |and d is the axis along which the pairs have zero spin projection, i.e. d · S pair ∆ = 0. This corresponds P For materials in which spin-orbit effects are strong, i.e. actinide and rare earth heavy fermion metals, the labels "spin-singlet" and "spin-triplet" do not refer simply to the eigenvalues of the spin operator for electrons, but rather a "pseudo-spin". The Kramers' degeneracy in zero-field guarantees that each p state is two-fold degenerate. These two states may be labeled by a pseudo-spin quantum number α. Thus, for many of our considerations the distinction between spin and pseudo-spin is not important [49,50,51].
to a spin-state that is an equal amplitude superposition of | ↑↑ and | ↓↓ spin states in the plane perpendicular to d. In bulk 3 He-A, which is stable over a narrow temperature range at high pressures, both the orbital and spin quantization axes are broken symmetry directions, aligned only by external walls or weak symmetry breaking fields. In particular, the very weak nuclear dipolar energy is minimized by orienting the nuclear spins in the orbital plane, i.e. d|| ± ℓ [53]. In thin films of superfluid 3 He scattering by the surface and substrate leads to substantial pair-breaking unless the pairs are confined to the x-y plane of the film [54]. This condition favors the A-phase with d||ℓ||ẑ over the B-phase for the entire pressure range (c.f. Fig. 3 in Ref. [21]). However, in this reduced geometry the planar phase, which is also an equal-spin pairing state, but with d z = 0 is degenerate with the A-phase in weak-coupling theory. It is unknown whether or not the A-phase or the planar phase is the stable phase in thin films, particularly at low pressures where strong-coupling effects are extremely small [55]. For the planar state the two spin states are aligned out of the plane of film and correlated with the orbital pairing state to form a state with zero total angular momentum projection alongẑ,

Broken Time-Inversion and Chirality
In contrast to the planar phase the A-phase spontaneously breaks both time-reversal symmetry and reflection symmetry in a plane normal to the film. This combination of broken 2D reflection symmetry and time-reversal symmetry is referred to as broken chiral symmetry. The structure of vortices in chiral superfluid 3 He-A are fundamentally different than vortices in non-chiral condensates such as the planar phase of superfluid 3 He. There are also candidates for triplet pairing in the class of strongly correlated electronic superconductors. For example, UPt 3 is believed to have an odd-parity order parameter belonging to the E 2u representation of D 6h with a low-temperature phase that is an equal-spin pairing state with d||ẑ, and breaks parity and time-reversal symmetry, i.e. ∆ = dp z (p x + ip y ) 2 . + The evidence and analysis favoring this identification comes from observations of a multiple superconducting phases, low-temperature transport properties and observation of anisotropic Pauli limiting, and is summarized in Refs. [56,22,15]. More recent small-angle neutron scattering studies of the flux lattice transitions also support this identification [57,58].
The layered superconductor, Sr 2 RuO 4 is also candidate for chiral spin-triplet pairing with a proposed order parameter that is essentially identical to that for a thin film of 3 He-A [59]. Strong evidence for unconventional pairing in Sr 2 RuO 4 is provided by + The higher symmetry of the hexagonal rotation group allows for a second two-dimensional E representation formed from the basis functions of the two B representations of D 4h . The odd-parity version, E 2u , with d||ẑ has basis functions dp z (p x ± ip y ) 2 .
the suppression of T c below 1.5 K by non-magnetic impurities [60], while the absence of a Hebel-Slichter peak in 1/T 1 T probed by NQR and Knight shift measurements [61] are cited in support a spin-triplet order parameter with the spin-quantization axis d||ẑ. Phase-sensitive measurements based on the current-phase relation for Josephson junctions coupling a conventional s-wave superconductor, Au-In alloy, and Sr 2 RuO 4 support the interpretation that Sr 2 RuO 4 is an odd-parity superconductor [62]. Studies of the excess current obtained from point contact spectroscopy agree quantitatively with a p-wave triplet pairing state [63]. Evidence for broken time-reversal symmetry of the superconducting order parameter is provided by µSR [64] and magneto-optical Kerr effect [65] experiments indicating that spontaneous currents and magnetic fields develop below T c in the Meissner state. However, the identification of broken chirality is controversial because local magnetic probes have not observed the fields associated with the domain-wall currents that are expected from broken chiral symmetry [66].

Vortex structure in chiral p-wave, spin-triplet condensates
Superconductors with a multi-component order parameter exhibit a variety of different vortex structures not possible in condensates belonging to a one-dimensional representation. These structures can lead to multiple superconducting phases as well as novel electronic and magnetic properties [17]. Theoretical studies of these vortex structures can provide important information for future experimental studies of vortex phases in triplet superconductors.
Theoretical analysis of the stability of the topologically stable vortex structures requires calculations of the free energy of these structures, typically taking into account the interactions with other vortices in the vortex lattice structure. Here we restrict our analysis and discussion to the internal structure of individual vortices appropriate to the low-field regime H ≃ H c 1 .
For any inhomogeneous superconducting state, and specifically for a vortex in a p-wave, spin-triplet superconductor or superfluid, the order parameter can be expanded in basis functions of the relevant irreducible representation, Γ, of the symmetry group, For chiral p-wave states appropriate for thin films of 3 He-A and Sr 2 RuO 4 with d|| z, the representation is two-dimensional and the basis functions, η Γ, ν (p f ) = d η ν (p f ) with η ± (p f ) =p x ±ip y , correspond to orbital angular momentum eigenstates with L orb z = ± . For 3 He-A these basis functions are exact to the extent that higher-order representations, e.g. f-waves with L z = ±1, can be negelected [67]. In the case of Sr 2 RuO 4 higher-order harmonics belonging to the same irreducible representation, E u , are possible. Typically, we can neglect these higher-order terms as the basic features of the vortex structure are obtained with the lowest order odd-parity harmonics.
Because time-reversal symmetry is broken for chiral, p-wave superconductors, there are two possible degenerate ground states defined by their spontaneous orbital motion.
The homogeneous ground state is spontaneously chosen between the two degenerate time-reversed states with L orb z = ± [46]. In the following sections we consider vortices in a single domain with internal orbital angular momentum L orb z = + . Vortices in the time-reversed ground state can be obtained simply by applying the time-inversion operation: ∆ → T * ∆.

Global and local winding numbers for vortices
Consider an isolated vortex with circulation m × 2π far from the vortex core, with m an integer. In this limit the amplitude approaches its bulk value but carries the global phase winding, Asymptotically the vortex is cylindrically symmetric, and is thus an eigenfunction of the generator, L z , for rotations about theẑ axis, * where L z = L cm z + L orb z is the sum of the z-component of the angular momentum operator for the center of mass (CM) of the pairs, L cm z = i ∂/∂φ with R = |R|(cos φ, sin φ), and the orbital angular momentum of the pairs relative to the CM, L cm z = i ∂/∂ϕp, where ϕp is the azimuthal angle of the relative momentum of the pairs, p f = p f (cos ϕp, cos ϕp). Thus, we have l = m + 1. At finite, but large distances from the vortex core, i.e. |R| ≫ ξ 0 , gradients of the asymptotic vortex order parameter induce the degenerate, time-reversed pairing state leading to an order parameter of the more general form, It is clear that the time-reversed phase with winding number p is nucleated in the vicinity of the vortex core. At large distances the vortex retains axial symmetry. This enforces a condition on the local winding number, p, of the time-reversed phase in the core, i.e. p − 1 = m + 1 [22]. This condition is exact for vortices in a pure p-wave chiral ground state such as 3 He-A. For chiral ground states in superconductors with discrete symmetry the condition on the local phase winding can be generalized. For example, for the discrete point group D 4h appropriate to Sr 2 RuO 4 the condition on the local phase winding becomes, where n = 0, ±1, ±2, . . . includes the effect of higher order harmonics comprising the ground state that are invariant under 4-fold rotations. If deviations from cylindrical symmetry are small we can neglect corrections coming from n = 0. For simplicity we Table 2. Table of vortex states for the chiral ground state with L orb z = + . Shown are the dominant order parameter, ∆ + , its global phase winding m, and the corresponding vortex core order parameter, ∆ − , and its local phase winding p. The color code describes the phase winding -a change in color corresponds to a change of the sign in the real part of the corresponding order parameter amplitude. The order parameter with m = 0 and p = +2 corresponds to a defect with no global phase winding, but with local phase winding near the defect. restrict our quantitative analysis to vortices with n = 0, but comment on the qualitative effects of the higher order harmonics where appropriate. Table 3.1 summarizes the lowest energy vortex states and their corresponding winding numbers. The figures show calculations of the dominant order parameter, ∆ + , as well as the sub-dominant order parameter, ∆ − , that develops in the core of the vortex. The color coding indicates the corresponding phase winding -a change in color corresponds to a change of the sign in the real part of the corresponding order parameter amplitude. The size of the circular symbols indicates the relative magnitude of the order parameter. We discuss each of the vortex states in more detail below. Figures 2, 3, and 6 provide additional quantitative results for the order parameter structure, the current distribution and local density of states for the vortices listed in Table 3.1. The magnitudes of the order parameter components, ∆ + and ∆ − , are shown in column 3 as a function of distance along the x axis through the vortex center. The local density of states, N (p f , R; ǫ) is shown for a family of trajectories parallel to the x axis with impact parameters y separated by ∆y = 0.78ξ 0 . Measurements of the local density of states using STM spectroscopy, which has been successfully applied to vortex studies in layered superconductors [68] including the high-T c cuprates [69], could provide valuable information on the vortex structure and pairing symmetry of Sr 2 RuO 4 .
Our calculations were performed with impurity scattering included in Born limit for a mean free path of ℓ = 10ξ 0 , where ξ 0 = v f /2πk B T c is the coherence length. The Fermi surface parameters are assumed to be isotropic, and the temperature was chosen to be T = 0.2 T c . For simplicity we also assumed the high-κ limit, where the penetration depth is large compared to the coherence length. Our calculations of the vortex structure are appropriate to the low-field limit near H c1 where vortices are well separated from each other, or for H < H c1 in the case of the defect with m = 0. The order parameter, impurity self energy, and the equilibrium spectra and current densities were obtained self consistently using the Riccati formulation of the quasiclassical transport theory with impurity and pairing self-energies.

Two types of singly quantized vortices
Broken time-reversal symmetry of the ground state leads to an interesting effect: vortices with opposite global phase winding numbers are inequivalent [17]. The ground state pairs have internal orbital angular momentum, L orb z , which can be either parallel or antiparallel to the external field H = Hẑ that nucleates a vortex and fixes the sign of the global phase winding. For H||L z = + ẑ the m = 1, p = 3 vortex is realized, whereas for the antiparallel case, −H||L z = + ẑ or H||L z = − ẑ, the m = −1, p = +1 vortex is nucleated. The structures of both vortices are shown in Fig. 2. Both vortices exhibit the large zero-energy Andreev bound state at the vortex center that is characteristic of odd winding number [11]. Note however the difference in magnitude of the induced components for m = ±1. The relative orientation of the field, H, and the spontaneously chosen direction for the internal orbital angular momentum of the pairs in the ground state, L orb z = ± ẑ, determines which type of singly quantized vortex is nucleated from the Meissner state at the lower critical field. The m = +1 vortex has a local phase winding of p = +3, and thus different core energy than the simpler m = −1 vortex with p = +1. As originally noted by Tokuyasu et al. [17], the difference in free energy for these two vortices leads to a splitting of the lower critical field, H c 1+ = H c 1− , for fields parallel (+) or antiparallel (-) to the internal orbital angular momentum of the pairs, L orb z [17]. Observation of this asymmetry would provide a direct signature of broken chirality in the ground state.
Ginzburg-Landau (GL) calculations by Tokoyasu et al. [17] also predicted the possibility of spontaneously broken axial symmetry by the cores of singly quantized vortices with m = ±1. In the GL theory whether or not the axial symmetry is spontaneously broken is determined by competition between the suppression of the condensation energy by the global phase winding and the internal Josephson phaselocking energy between the two time-reversed order parameter amplitudes. If the GL coefficient, β 2 , that determines the coupling energy, F c = 4β 2 |∆ + | 2 |∆ − | 2 , is sufficiently weak compared with the condensation terms, β 1 [|∆ + | 4 + |∆ − | 4 ], i.e. 0 < β 2 ≤ 1 4 β 1 then axial symmetry is spontaneously broken by the vortex core. Our results for singly quantized vortices, which do not break axial symmetry at any temperature, also agree with those of Ref. [17] since our theory reduces to the GL equations with the weakcoupling value of β 2 /β 1 = 1 2 in the limit T ≈ T c . Other authors have also considered differences between singly quantized vortices with opposite phase windings relative to the chirality of the ground state. The authors of Ref. [70] report results for the field dependence of the vortex structure in chiral p-wave superconductors. They argue that H c 2 differs for states with opposite chirality, i.e. L orb z = ± . However, this is incorrect. The upper critical field represents a second-order transition between the normal state and the superconducting vortex state. Since chiral symmetry is unbroken in the normal state, the second-order instability field, H c 2 , is necessarily the same for H|| +ẑ or H|| −ẑ. Furthermore, at H c 2 both chiral components of the 2D representation are nucleated. As the field is reduced, the vortex density decreases and eventually a single chiral domain may remain, at least in the idealized limit of a perfect crystal with no vortex pinning. The chirality of the resulting ground state in the limit H → 0 will be determined by the direction of H, but H c 2 is the same for either field orientation. The situation is completely different for H c 1 . In this case a vortex enters the superconductor from the Meissner state which spontaneously breaks chiral symmetry. The relative sign of the phase winding of the vortex and the chirality of the pre-existing Meissner phase lead to different vortex core states nucleating depending on whether H|| +ẑ or H|| −ẑ, and a corresponding asymmetry in H c 1 .
Recently, Yokoyama et al. [71] proposed that broken chiral symmetry should be observable as an asymmetry in the surface density of states (SDOS) at the Fermi level when a vortex is situated within a coherence length of the surface of a chiral p-wave superconductor. These authors predict a large anomaly in the SDOS depending on the sign of the chirality based on a non-self-consistent model for the vortex core which omits the internal structure of the vortex cores of chiral superconductors described above.
Similarly, the authors of Ref. [72] report a dramatic difference in the broadening of the energy levels in the vortex core for the m = +1 and m = −1 vortices.♯ Their analysis is based on an assumed form for the vortex amplitude that does not include the induced component of the order parameter in the core with time-reversed chirality. Our results, which are based on fully self-consistent solutions for the order parameter, single-particle self-energy and spectral density, show only a minor reduction in the spectral width of the m = −1 core state compared to that for the m = +1 core state. This is consistent with earlier calculations by Hess et al. [73,74], who found a small asymmetry in the dispersion of the Andreev bound states for self-consistent solutions for singly quantized vortices in the chiral ground state.
The discrete symmetry of the point group is not expected to significantly change the structure of the m = −1, p = 1 vortex since the next relevant higher order components are m = −1, p = −3, +5 . . .. However, the m = +1, p = +3 vortex allows a higher-order subdominant component with p = −1, which may have a significant amplitude since it will be localized relatively close to the vortex center. The corrections to the vortex core resulting from higher order harmonics to the basis functions allowed by the discrete point group symmetry will reduce the difference in free energy between the two types of singly quantized vortices. Takigawa et al. [75] have considered the extreme anisotropic limit for the structure of singly quantized vortices and the local DOS by solving the ♯ The m = +1 vortex in the chiral ground state with L orb z = + ẑ corresponds to the L z = +2 vortex of Ref. [72], while the m = −1 vortex in the same chiral ground state corresponds to their L z = 0 vortex.  Table 3.1 and the size of the symbols reflects the magnitude of the order parameter. In the third column we show both order parameter components, ∆ + (black) and ∆ − (red), as a function of distance along a trajectory through the vortex center. In the fourth column the local density of states, N (x = 0, y; ǫ), is shown as a function of energy and "impact parameter", y, from the the vortex center for y = −12.5ξ 0 to y = 12.5ξ 0 and a spacing of ∆y = 0.78 ξ 0 . The red spectrum showing the zero-energy core state is evaluated at x = 0, y = 0.
Bogoliubov equations numerically on a lattice with nearest neighbor (NN) hopping as well as local NN pairing interactions.

Doubly quantized vortices
As in the case of singly quantized vortices, there are two in-equivalent doubly quantized vortices: (i) m = 2, p = 4 and (ii) m = −2, p = 0 as shown in Fig. 3. The latter case describes a doubly quantized vortex with a nearly homogeneous core amplitude, i.e. a "coreless vortex". Since there is no local phase winding associated with the subdominant amplitude, it develops to the asymptotic value for the dominant order parameter component as indicated in the third column of Fig. 3. This is similar to a domain wall separating the two degenerate order parameter components, except that the subdominant component is restricted to the core (see Fig. 4). The large subdominant component which "fills in" the core of the doubly quantized vortex substantially lowers the free energy for the m = −2, p = 0 vortex by recovering nearly all the lost condensation energy from the core of the dominant component. This low core energy is reflected in the local density of states shown in the lower panel of column 4 in Fig.  3; there are very few low-energy, sub-gap excitations associated with pair-breaking in the vortex core. The large reduction in the core energy of the doubly quantized vortex compared to that of the singly quantized vortex allows for stable lattices of  doubly quantized vortices at sufficiently high magnetic fields. The possibility of stable lattices of doubly quantized vortices in superconductors with broken chiral symmetry was investigated within GL theory. It was shown that doubly quantized vortices can be stabilized compared to the singly quantized vortices over a wide range of fields below H c2 [76]. The axial symmetry of the current distribution (shown in the lower panel of Fig. 4) leads to isotropic interactions between these doubly quantized vortices, and thus a hexagonal lattice structure.
Broken Axial Symmetry in the Vortex Core The calculation of the structure of the doubly quantized vortex with m = +2 and p = +4 shown in Fig. 3 is for a relatively small area grid. The result shows an axially symmetric structure in which the local phase winding of p = +4 is concentrated as a multiply-quantized axially symmetric vortex. However, this structure is unstable. The stable solution for this vortex spontaneously breaks axial symmetry, with the p = +4 vortex in the time-reversed phase dissociating into four singly quantized vortices arranged around the edge of the time-reversed core amplitude as shown in the top-right panel of Fig. 5. The breaking of C ∞ to C 4 symmetry is energetically preferred because dissociation of the p = +4 vortex allows the timereversed amplitude in the core to recover to the bulk value and restore most of the lost condensation energy of the dominant phase. This reduces the energy splitting between the two types of doubly quantized vortices and implies that the region of stability of doubly quantized vortices will be comparable for either field orientation, H||L z = + ẑ or −H||L z = + ẑ. However, the axial anisotropy of the current distribution of the doubly quantized vortex (lower panel of Fig. 5) is expected to lead to a square vortex lattice of doubly quantized vortices for H||L z = + ẑ, compared to a triangular lattice for −H||L z = + ẑ. Finally, let us consider the effects on vortex structure from additional order parameter components which are allowed by discrete point symmetry of Fermi surface. For the axially symmetric double quantum vortex with m = −2 and p = 0 the next order harmonics are those with p = ±4. The magnitudes will be very small because both amplitudes, ∆ + and ∆ − , are large and there is no energetic advantage to inducing these components. Thus, we can safely neglect these higher harmonics. The same reasoning applies for the vortex with m = +2 and p = +4. Both the dominant and core amplitudes are large and suppress the higher harmonics, even the higher harmonic with n = −4 (i.e. p = 0). This would not be the case if the axial symmetry were not spontaneously broken.

Structure of inhomogeneities in the Meissner phase
In addition to novel vortex states, broken chirality leads to defect structures that carry current, but have no global phase winding. The example of a cylindrical normal metallic inclusion is included in Table 3.1. This defect has no global phase winding, i.e. m = 0, but necessarily has a local phase winding of p = +2 for the induced order parameter. Thus, in contrast to vortices that are stabilized by a field, this defect can exist in the Meissner phase. Such defects can be engineered, or as in the case of Sr 2 RuO 4 , Ru ions may precipitate and form metallic inclusions embedded in the superconducting host material.
The ground state of a chiral, p-wave superconductor will be either thep x + ip y state, or thep x − ip y state, although the actual zero-field superconducting phase may include a number of domains of different time-reversed states separated by domain walls. Here we consider a normal metallic inclusion of radius r incl = 0.78ξ 0 embedded in the homogeneousp x + ip y ground state, or at least several coherence lengths away from a domain wall. We assume a simple model for a metallic inclusion in which electrons and holes are perfectly transmitted across the (NS) interface between the inclusion and the superconducting material when the latter is in its normal state. We also assume there is no pairing interaction within the metallic inclusion. Thus, a mean field order parameter does not develop inside the inclusion, however pairing correlations and associated spontaneous currents do develop within the inclusion as a result of the proximity effect and the high transmission of the NS interface. More detailed models, including specific models of Ru inclusions, which incorporate the finite reflectivity of the NS interface, as well as pairing attraction within the metallic inclusion can be implemented, but are outside the scope of this article.  Figure 6. The structure of a normal metallic inclusion embedded in the (p x + ip y ) phase. The inclusion is located in the center in column 1 and has a radius r incl = 0.78ξ 0 . The columns and notation are the same as that of Fig. 2. Note (i) the phase winding of the time-reversed order parameter (column 2) and (ii) the shallow bound states induced by the inhomogeneous order parameter near the boundary of the metallic inclusion (column 4). Fig. 6 is the time-reversed component of the order parameter, with a phase winding of 4π, that nucleates near the metallic inclusion. Broken timereversal symmetry is revealed by the appearance of spontaneous supercurrents that flow in and around the metallic inclusion as shown in the top-right panel of Fig. 7. These currents generate a local magnetic field, b = b(x, y)ẑ, which oscillates in sign as one moves radially away from the defect on the scale of the coherence length, ξ 0 , shown in top-left panel of Fig. 7. The spatial average of the field is zero, but the locally varying field is non-zero and should be observable, e.g. as broadening of the µSR linewidth below T c . We also note that the current and field induced by a nonmagnetic point impurity in a chiral p-wave superconductor [77] is similar to that of a mesoscopic metallic inclusion, however the origin and interpretation of the structure appears different. Complex current and field distributions may also be generated in unconventional, but non-chiral, superconductors by magnetic and spin-orbit scattering impurities [19].

Shown in column 2 of
Finally we note that the leading higher order harmonic corrections resulting from discrete C 4 lattice symmetry are expected to be weak for a small cylindrically symmetric defect. The most important higher harmonic will be a subdominant component with p = −2, but it is weaker than the leading p-wave harmonic with p = +2. All other harmonics can be neglected for small inclusions because the higher winding numbers will force the induced amplitude to be further from the vortex core and thus more effectively suppressed by the dominant m = 0 order parameter. A possible exception may occur for large inclusions, r incl ≫ ξ 0 . In this limit the harmonic with p = 6 may be more important than the p = −2 component.

Conclusions
In conclusion we have analyzed and carried out self-consistent calculations of the equilibrium structure and excitation spectrum of vortices in layered, spin-triplet, p-wave superconductors and superfluid 3 He films with spontaneously broken chiral symmetry in the ground state. We show that the cores of vortices contain the time-reversed phase with a local phase winding that differs by ±2 × 2π depending on the chirality, i.e. L orb z = ± , of the ground state. This implies that vortices with equal but opposite phase winding are inequivalent, and thus the lower critical field for nucleation of vortices in layered, chiral p-wave superconductors depends on the relative orientation of the field and the orbital angular momentum of the ground-state pairs. Axially symmetric doubly quantized vortices with zero phase winding in the core are predicted. A lattice of these vortices is energetically favored at sufficiently high fields for −H||L z = ẑ. For the opposite phase winding, i.e. H||L z = +ẑ, the doubly quantized vortex contains four circulation quanta for the core amplitude. This vortex spontaneously breaks axial symmetry by dissociation of the core vorticity into four singly quantized vortices in the time-reversed order parameter induced in the core. A lattice of these vortices is also expected to be energetically favored at sufficiently high fields with a square lattice structure. In addition to vortex states, inhomogeneities in the distribution of non-magnetic impurities, or metallic inclusions embedded in the ground state lead to spontaneous supercurrents and magnetic fields induced in the Meissner phase, which are localized near the inhomogeneity. These currents are carried by bound electronic states associated with an induced order parameter with phase winding of 4π. The magnetic field and current distribution vary on the scale of the coherence length and are expected to be observable by sufficiently small local magnetic probes, or as a broadening of the µSR linewidth below T c .