Defect bulk-boundary correspondence of topological skyrmion phases of matter

Unpaired Majorana zero-modes are central to topological quantum computation schemes as building blocks of topological qubits, and are therefore under intense experimental and theoretical investigation. Their generalizations to parafermions and Fibonacci anyons are also of great interest, in particular for universal quantum computation schemes. In this work, we find a different generalization of Majorana zero-modes in effectively non-interacting systems, which are zero-energy bound states that exhibit a cross structure -- two straight, perpendicular lines in the complex plane -- composed of the complex number entries of the zero-mode wavefunction on a lattice, rather than a single straight line formed by complex number entries of the wavefunction on a lattice as in the case of an unpaired Majorana zero-mode. These cross zero-modes are realized for topological skyrmion phases under certain open boundary conditions when their characteristic momentum-space spin textures trap topological defects. They therefore serve as a second type of bulk-boundary correspondence for the topological skyrmion phases. In the process of characterizing this defect bulk-boundary correspondence, we develop recipes for constructing physically-relevant model Hamiltonians for topological skyrmion phases, efficient methods for computing the skyrmion number, and introduce three-dimensional topological skyrmion phases into the literature.


I Introduction
Discovery of the integer and fractional quantum Hall effects in the early 1980's pushed condensed matter physics beyond the Landau paradigm 1 of symmetry-breaking phases of matter, towards additional characterization of matter based on its topology. Topological phases of matter, or those phases of matter with some observables linked to topological invariants and therefore unaffected by sufficiently small perturbations, have profound consequences. They realize quasiparticles useful for exceptionally robust quantum computation schemes 2 , and their topology can manifest through bulk-boundary correspondence: a non-trivial value for a topological invariant of a material bulk yields topologically-protected, gapless boundary states 3-5 very promising for applications such as spintronics 6 and topologically-protected quantum computation 7,8 .
Study of topological phases has greatly expanded since the discovery of the quantum spin Hall insulator, a topological phase of matter realized without an applied magnetic field or magnetic order, presenting the possibility of far broader applications than topological phases studied earlier 9 . Efforts began in earnest to search for other previously-unidentified topological phases of matter, leading to the three-dimensional topological insulator protected by time-reversal symmetry 10 , and later the ten-fold way classification scheme for topological phases of matter 11,12 . Much more recently, efforts to expand the ten-fold way to include bosonic topological phases 13 , crystalline topological phases [14][15][16][17][18] , the broader set of short-range and long-range entangled topological phases 19,20 , topology of various quasiparticles 21? -25 and driven [26][27][28][29][30] and open [31][32][33][34] systems reveal the foundational role of topology in condensed matter physics.
In all of these topological phases of matter, however, topology is associated with the full set of degrees of freedom of a system. This reflects the characterization of many topological phases based on projectors onto occupied states 11,12 . Recent work has shown, however, that topological phases are also realized for subsets of the degrees of freedom of the occupied states. Specifically, past work showed non-trivial topology is also associated with just the spin degree of freedom (DOF) of the ground-state in symmetry-protected topological phases even when the ground-state is characterized by a larger set of degrees of freedom 35 , independent of the previously-known topology associated with the full set of degrees of freedom of the ground-state.
These topological phases of the spin degree of freedom, known as topological skyrmion phases of matter, serve as momentum-space duals of the widely-studied skyrmion magnetic orders in real-space [36][37][38][39] . As these skyrmions are realized in momentum-space and their symmetry-protection is robust against crystal-field splitting, however, they are topological in a stricter sense, exhibiting a bulk-boundary correspondence upon performing a partial trace over all degrees of freedom of a system except spin. As these topological phases are protected by a generalized particle-hole symmetry present in centrosymmetric superconductors, such signatures may be important for experiments and understanding of exotic mechanisms of superconductivity. The consequences of this nontrivial topology are even more significant, however, yielding novel topological phase transitions even in effectively noninteracting systems while the symmetry protecting this spin topology is preserved, as the topology is protected by a spin magnitude gap rather than an energy gap. Thus, given that the topological skyrmion phases are experimentally-relevant and generalize foundational concepts of topological condensed matter physics, better understanding of these phases will yield valuable insights.
In this work, we examine these topological skyrmion phases of matter, revealing they exhibit a second kind of bulkboundary correspondence due to the momentum-space spin textures trapping defects, while also introducing the threedimensional topological skyrmion phases of matter into the literature. These results agree with past work on topological classification of defects by Teo and Kane 40 , but there are surprising additional consequences of topological defects in momentum-space. Past work shows unpaired Majorana zero-modes may be realized for certain boundary conditions when spin skyrmions are realized in momentum-space for the Hopf insulator [41][42][43][44] . We extend these results to the topological skyrmion phases of matter, showing these momentum-space spin skyrmions realize more exotic topologically-protected zero-modes. In identifying this second defect-based bulkboundary correspondence of the topological skyrmion phases of matter, we provide strong evidence that the consequences of this physics, including the type-II topological phase transition of the topological skyrmion phases, must be explored in detail to better understand the foundations of topological condensed matter physics.
The overview of the manuscript is as follows: in Sec. II, we review the topological skyrmion phases and then expand on earlier results by presenting methods for constructing toy models-including those corresponding to Bogoliubov de Gennes Hamiltonians for superconductors-and highlyefficient calculation of the topological invariant. We use these techniques to characterize the models numerically, studying their phase diagrams and momentum-space spin textures. In Sec. III, we then extend our discussion by introducing threedimensional topological skyrmion phases using the construction methods of Sec. II. Using these methods, we construct toy models for the three-dimensional skyrmion phase from the Hopf insulator 43,44 . We then show these three-dimensional spin skyrmions trap momentum-space defects, which realize remarkable generalizations of unpaired Majorana zero-modes for certain boundary conditions. We characterize these socalled cross zero-modes numerically and analytically, showing how they relate to the unpaired Majorana zero-mode. Finally, in Sec. IV, we provide some summarizing remarks and discuss possible extensions for future research.

A Section overview
In this section, we first provide a brief overview of the topological skyrmion phases and the standard method for computing the topological invariant characterizing them. We then introduce the simplest, or "atomic", symmetry-allowed toy model Bloch Hamiltonians that can realize this non-trivial topology. Based on these atomic toy models, we then present some general design principles for constructing toy model Bloch Hamiltonians with four bands that can realize the chiral phase. We then introduce a far more efficient method for computing the topological invariant and discuss situations where it is most useful. Finally, we compute phase diagrams to characterize the topology of some canonical four-band toy models exhibiting the chiral topological skyrmion phase.

B Review of the chiral topological skyrmion phase
The chiral topological skyrmion phase is a topological phase of matter introduced in past work 35 , which is characterized by formation of a baby skyrmion in the ground state spin expectation value texture over the two-dimensional Brillouin zone, even when systems have more degrees of freedom than just spin. The phase is realized in Bloch Hamiltonians with an even number of bands greater than two, and has previously been studied in a tight-binding model relevant to superconducting Sr 2 RuO 4 . This non-trivial topology results when the Bloch Hamiltonian is invariant under a generalized particlehole transformation corresponding to operator C ′ , which is the product of the particle-hole operator C (with C 2 = −1), and the spatial inversion operator, I. This symmetry ensures that the Hamiltonian serves as a mapping from the Brillouin zone to a space acted on by the quotient Sp(2N )/U(N ), where Sp(2N ) is the compact symplectic group with Sp(2) (N = 1) isomorphic to the special unitary group SU(2). One physical quantity acted upon by such a quotient is the spin angular momentum. As a result, systems with C ′ symmetry can realize topological phases characterized by the skyrmion number Q computed as an integral over the Brillouin zone of the normalized ground state spin expectation valueŜ(k) = S(k)/|S(k)|, where S(k) = n∈filled ⟨u n (k)|S|u n (k)⟩ and the spin operator S reads is a set of Pauli matrices acting in generalized particle-hole (spin-1 2 ) space. τ 0 and σ 0 are 2 × 2 identity matrices.
We may express the skyrmion number as The integrand here is one component of an object we call the skyrmion curvature Ω Q , expressed as where the k th component is defined as where ϵ ijk is the Levi-Civita tensor and Einstein summation is implicit.

C
Methods for constructing toy model Hamiltonians realizing the 2D chiral topological skyrmion phase In order to construct toy models realizing the chiral topological skyrmion phase, we require that they serve as mappings from the Brillouin zone to the space of ground state spin expectation values corresponding to non-trivial homotopy group π 2 (Sp(2N )/U(N )) = Z. As discussed in earlier work, this requires the Hamiltonian be invariant under a C ′ transformation 45 . For a 2N -band model, the C ′ symmetry operates on fermion operators as The Hamiltonian matrix representation must therefore satisfy The form of the Hamiltonian for N = 2 (corresponding to four bands) is therefore restricted to include some or all of the ten checked terms in TABLE I, expressed as Kronecker products of Pauli matrices τ i σ j = τ i ⊗ σ j : We further classify these ten terms based on how they affect the ground state spin expectation value for half-filling. To simplify the problem, we first consider each term individually, computing the spin expectation value for each τ i σ j term, of the ten terms allowed by C ′ symmetry, separately. We find that three of the ten terms yield non-zero values for the x-, yor z-components of the spin when considered alone: specifically, τ 3 σ 1 alone yields non-zero ⟨S x ⟩, τ 0 σ 2 alone yields non-zero ⟨S y ⟩, and τ 3 σ 3 alone yields non-zero ⟨S z ⟩. That is, there exists one term τ α σ β with four eigenvectors labeled by i ∈ {1, 2, 3, 4} as {ψ i,αβ } satisfying the eigen problem τ α σ β ψ i,αβ = E i,αβ ψ i,αβ and eigenvalues {E i,αβ } satisfying E 1,αβ ≤ E 2,αβ ≤ E 3,αβ ≤ E 4,αβ , such that ⟨S γ ⟩ = i∈{1,2} ψ † i,αβ S γ ψ † i,αβ ̸ = 0 for each γ ∈ {x, y, z}. As the chiral skyrmion phase requires realization of a skyrmionic spin texture in the Brillouin zone and therefore finite x-, yand z-components of the ground state spin texture somewhere in the Brillouin zone, the simplest toy model for a non-trivial chiral skyrmion phase must therefore include these three terms. This furthermore suggests that the simplest way to construct various four-band Bloch Hamiltonians restricted to these three terms is as where {d i (k)} are the momentum-dependent functions of a two-band model for a Chern insulator expressed as The other seven symmetry-allowed terms may be classified depending on how they alter the x-, y-, or z-components of the ground state spin texture. Specifically, a spin component can be finite when one of these seven terms is multiplied by the same constant as one of the three essential terms, or zero. We can therefore construct TABLE II, marking each of the additional seven terms as yielding F (finite) or Z (zero) spincomponent when it is multiplied by the same prefactor as one of the three essential terms. This forms the basis for construct-

D Relation between the Chern number and the Skyrmion number in a special model
Here we investigate the relation between the Chern number and the skyrmion number in a special case.
The Hall conductivity in an ideal lattice evaluated directly from the Kubo formula reads where m, n are band indices, ω and η are set to zero in the dc clean limit, the velocity operator v a = ∂H(k)/ℏ∂k a (a = x, y) in Bloch basis, and V is the volume of the crystal in the real space. For a 4-band model described by H(k) = α d α (k)Γ α , we can express the velocity operator explicitly as And it follows that the matrix element for the velocity operator reads where ⟨Γ α (k)⟩ nm ≡ ⟨u n (k)|Γ α |u m (k)⟩. Using the explicit form of velocity matrix elements, the dc Hall conductivity with C the well-known Chern number. We note that the second equation is valid only when the chemical potential is inside the gap.
(2), we may also compute it using a far more efficient method based on projectors, developed for computing the Chern number in past work 46,47 . The idea of this earlier work was that computing a Berry phase directly from summation of Berry curvature is more costly than computing it as a sum of infinitesimal Berry phases, each computed using a discretized form of the Berry connection. These past works found that the latter approach converged to the correct, integer-quantized values for Chern numbers using a far coarser mesh of k-points in the Brillouin zone than did the former. In analogy to this approach, we want to compute the skyrmion number with a connection, rather than the skyrmion curvature defined previously. We therefore first present an efficient method for computing this connection starting from the ground-state spin expectation value texture over the Brillouin zone.
To compute this skyrmion connection, we first construct an eigenstate of effective two-band model from the numerically computed ground state spin expectation value ⟨Ŝ(k)⟩ of our four-band model.
Assuming half-filling of this effective two-band model, we may then write the form of the occupied eigenstate |ψ − (k)⟩ analytically in terms of {d S,i } as |ψ − (k)⟩ may then be used to compute the skyrmion number as the Chern number for the lower band C − of H eff , given that the skyrmion number and Chern number are equal in twoband models at half-filling.
We evaluate C − using a discrete lattice corresponding to the discrete values of the crystal momentum in a finite sample with periodic boundary conditions. The Chern number of the lower band may then be evaluated as where A n k,µ = ⟨ψ − (k)|ψ − (k + µ⟩ is the Berry connection field of the two-band model and the {µ} are the nearestneighbor vectors of the square momentum-space lattice. In the context of this effective two-band model, we call this the i th component of the skyrmion connection A Q , given that this effective two-band model is constructed from the ground state spin texture of a model for a chiral topological skyrmion insulator phase with more than two bands. This connection is more generally expressed, for the − band, as with the i th component expressed as The corresponding skyrmion curvature, which is the curl of this skyrmion connection, is the integrand of Eq. 3.
While this approach of relating the skyrmion number in a 2N -band model to the Chern number of a two-band model to make use of additional computation methods is generally not needed at the four-band level for skyrmion textures with topological charge sufficiently small in magnitude, it becomes invaluable for higher-magnitude topological charges or skyrmion numbers in models with more bands. As we show later, this effective connection may be used to compute the three-dimensional skyrmion number as well.

F
Characterization of canonical toy models Here, we construct some example toy models for N = 2 chiral topological skyrmion phases in two dimensions. We first present some example Hamiltonians to illustrate the decoupling of the total Chern number and skyrmion number at the four-band level by computing phase diagrams. We then discuss specific signatures of non-trivial skyrmion number at the four-band level beyond the momentum-space skyrmionic spin texture itself: a contribution to the Hall conductivity that may be expressed in terms of a quantized skyrmion number in C ′invariant models, and finite spin current signatures required for a system to possess non-trivial skyrmion number.
We first consider a model with quantized total Chern number and quantized skyrmion number, corresponding to co-existing stable Chern insulator and chiral topological skyrmion insulator phases, by restricting our Hamiltonian to only the three special terms in TABLE II. We then show that a model constructed using three of the seven other terms in the table may realize a Chern insulator phase characterized by a stable non-trivial total Chern number at half-filling but be topologically unstable according to the skyrmion number. After this, we consider toy models for superconductors written in Bogoliubov de Gennes (BdG) formalism, discussing constraints on the Hamiltonian required for it to describe topological skyrmion phases of matter in this physical setting and consider a representative case with spin singlet pairing and another case with spin-triplet pairing.

1
Toy model with coexisting Chern insulator and skyrmion insulator phases As the three special terms mentioned above can each support a non-zero spin component, we first consider the following Bloch Hamiltonian: To characterize the topology of this model, we compute the total Chern number C and skyrmion number Q assuming half filling of the bands, as well as the minimum direct gap between the bands second-lowest and third-lowest in energy over the Brillouin zone, d 2−1 , and the minimum magnitude of the ground state spin expectation value over the Brillouin zone |⟨S⟩| min , to obtain the phase diagrams shown in Fig. 1. The total Chern number C, shown in Fig. 1(a), and the skyrmion number Q, shown in Fig. 1(b), are each integer-quantized. We find a relationship between them of C = −2Q, in agreement  with non-trivial skyrmion number Q. Minimum direct gap between the bands second-lowest and third-lowest in energy over the first Brillouin zone is shown in (a). x, y, and z components of the normalized ground-state spin expectation value, ⟨Sx⟩, ⟨Sy⟩, and ⟨Sz⟩, over the first Brillouin zone are shown in (b), (c), and (d), respectively, and correspond to winding of the spin in agreement with the non-trivial skyrmion number value shown in Fig. 1 with Eqn. (15) for the atomic model. The minimum direct gap between the second-lowest and third-lowest bands in energy, d 2−1 , is shown in Fig. 1(c), and indicates that all topological phase transitions are type-I: they occur with a closing of a direct gap. Type-II topological phase transitions are not expected at the four-band level, due to restrictions on DOFs of the Bloch Hamiltonian other than the spin-1/2 DOF. We also find the minimum magnitude of the ground state spin expectation value over the Brillouin zone is finite everywhere in the phase diagram, as shown in Fig. 1(d), as required for well- defined and integer-quantized Q.
We examine a point in the phase diagram corresponding to β = κ = 1 in greater detail, computing the minimum direct gap between the second-lowest band in energy and the thirdlowest band in energy, d 2−1 , over the Brillouin zone, as well as the momentum-space spin texture, as shown in Fig. 2. d 2−1 over the Brillouin zone is shown in Fig. 2 (a) and confirmed to be finite everywhere the Brillouin zone. It also reflects the four-fold rotational symmetry of the model. Each of the three components of the spin expectation value are shown in Fig. 2(b), (c) and (d) respectively. They collectively show the winding of the spin vector in momentum space corresponding to a baby skyrmion.

2
Toy model with well-defined total Chern number and ill-defined skyrmion number To further demonstrate that the total Chern number and skyrmion number are decoupled at the four-band level, we consider a model constructed from a two-band Chern insulator Hamiltonian, which has zero magnitude for the ground state spin expectation value everywhere in the Brillouin zone. We choose the following Bloch Hamiltonian: and characterize it at half-filling by computing phase diagrams for the total Chern number C, skyrmion number Q, minimum direct gap between the bands second and third-lowest in energy, d 2−1 , and the minimum ground state spin magnitude over the Brillouin zone, |⟨S⟩| shown in Fig. 3. The phase diagram for the total Chern number of this Hamiltonian (Eq. (22)) in Fig. 3 (a) is identical to that of the total Chern number for the Hamiltonian Eq. (21). However, the skyrmion number is uniformly zero as shown in Fig. 3(b), because the x and y spin components do not wind, as shown in Fig. 3 (b) and (c), respectively. Additionally, the skyrmion number is actually ill-defined, as the minimum magnitude of the ground state spin expectation value over the Brillouin zone is zero as shown in Fig. 3(d). Specifically, the normalized xand y-components of the spin are ill-defined in the Brillouin zone as the unnormalized components are zero, corresponding to noise shown in Fig. 4(b) and (c). The z-component of the ground state spin expectation value also goes to zero at certain k points as shown in Fig. 4(d), resulting in zero magnitude for the spin expectation value at these momenta.

3
Toy models for chiral skyrmion superconductor phases As the generalized particle-hole operator C ′ is the product of the particle-hole operator C and the spatial inversion operator I, we naturally expect the skymion phases in centrosymmetric superconductors. Indeed, this has already been confirmed for tight-binding models of Sr 2 RuO 4 48 . To better understand when topological skyrmion phases are possible in such physical settings, we first consider a generic four-band model for a superconductor written in BdG formalism. The Hamiltonian H takes the form where and Here, c ks is the annihilation operator of electrons with momentum k = (k x , k y ) and spin angular momentum DOF s ∈ {↑, ↓}. E(k) is the Hamiltonian of the normal state, and ∆(k) is the gap function of the superconductor. This model possesses C ′ symmetry when −ϵ T (−k) = −ϵ * (k), which is satisfied for ϵ(k) even in momentum k. Therefore, to construct some representative example chiral topological skyrmion superconductor Hamiltonains, we take the normal state Hamiltonian to describe a Chern insulator that is even in momentum. One suitable example model is the following, from Sticlet et al. 49 : which realizes topologically-distinct insulating phases characterized by Chern numbers of 0, ±2 for the lower band. We may then combine this with a pairing matrix ∆(k) describing either spin singlet pairing ∆ s (k) or spin-triplet pairing ∆ t (k). We take with the d-vector of the spin-triplet pairing term here chosen specifically to be d(k) = sin(k y )x − sin(k x )ŷ as one example. This d-vector choice has previously been proposed as characterizing Sr 2 RuO 4 in the high-field phase 48 .
We may now compute phase diagrams for these two models characterizing spin-singlet and spin-triplet superconductivity, respectively. We first consider the case of spin-singlet pairing.
To characterize the topology of the model, we compute the total Chern number C and skyrmion number Q assuming two bands are filled, as well as the minimum direct gap between the third-lowest and second-lowest bands, d 2−1 , and the minimum magnitude of the ground state spin expectation value The total Chern number C phase diagram shown in Fig. 5 (a) may be divided into three distinct parts: a topologically stable region with total Chern number of 4 shown in red, which has an hourglass shape, two triangular regions where the Chern number is well-quantized at zero shown in green, and regions where the total Chern number is highly unstable. In comparison, the phase diagram for the skyrmion number in Fig. 5 (b) possesses a tan hourglass region with Q = 2 where the total Chern number C = 4. This is further evidence that the skyrmion number differs from the total Chern number by a factor of 2 at the four band level when they are each topologically stable. There are then two triangular regions with the skyrmion number very well-quantized at zero, and regions around the hourglass where the skyrmion number appears unquantized but close to zero.
To better understand the model, we also examine the phase diagram for the minimum direct gap d 2−1 , shown in Fig. 5 (c). We see that it is finite only for the hourglass region, explaining the instability observed in Fig. 5 (a) and (b). In addition, we compute the phase diagram for the minimum magnitude of the ground state spin expectation value, |⟨S⟩|. We see that it is finite in the hourglass region as required for a topologicallystable value of the skyrmion number, but zero in the region where the skyrmion number is not well-quantized. This is as expected, since the skyrmion number is computed with the normalized ground state spin expectation value, and is therefore unstable when the spin is not normalizable somewhere in the Brillouin zone.
We now consider the case of spin-triplet pairing. By tuning the parameters α and ∆ 0 , we again compute the total Chern number C and skyrmion number Q assuming half-filling of bands, as well as the minimum direct gap over the Brillouin zone d 2−1 and the minimum ground state spin expectation value over the Brillouin zone |⟨S⟩|.
The total Chern number phase diagram, shown in Fig. 6  (a), is topologically-stable almost everywhere in the phase diagram at values C = ±4. The skyrmion number phase diagram, shown in Fig. 6 (b), is similar in form, yet the tan regions with Q = 2 appear far more stable than the blue regions with Q ≈ −2. We see, again, that topological phase transitions occur where the minimum direct gap over the Brillouin zone between the two middle bands, d 2−1 , goes to zero, as shown in Fig. 6 (c). We gain insight into the difference between the Q = 2 regions and the Q = −2 regions by examining the phase diagram for the minimum ground state spin expectation value over the Brillouin zone, |⟨S⟩|: regions with Q = 2 have finite spin everywhere in the BZ, while |⟨S⟩| is zero in magnitude somewhere the BZ for Q = −2 regions. Thus these regions are topologically unstable according to the skyrmion number while being topologically stable according to the total Chern number, further indicating that these two topological invariants are decoupled from one another at the four band level. Remarkably, however, the ground state spin texture nearly forms skyrmions even in these cases where the spin is somewhere zero in magnitude in the BZ, similar to results for the mirror subsectors of the previously-considered tight-binding model for Sr 2 RuO 4 48 . Noticeably, only type-I topological phase transitions occur in each model, that is, only topological phase transitions that occur with a closing of the direct gap between the two middle bands. This results from the four-band C ′ -invariant Hamiltonians lacking a suitable DOF, such as orbital angular mo- mentum, into which spin angular momentum can flow. Although one could interpret the Bloch Hamiltonian as corresponding to a system with a spin half DOF and a two-fold orbital DOF, the C ′ invariance of the model restricts the orbital DOF such that only the type-I topological phase transition occurs at the four-band level. In Hamiltonians with six or more bands, however, topology of an orbital DOF is no longer ensured by C ′ symmetry as discussed in work introducing the topological skyrmion phases of matter and observed for tightbinding models describing superconducting Sr 2 RuO 4 48 , making the type-II topological phase transition possible.
We fix a particular choice of parameters ∆ 0 = 1, α = 2 to study the spin textures in the Q = −2 skyrmion phase. As shown in Fig. 7   In this second section of the paper, we introduce the threedimensional counterpart of the two-dimensional chiral topological skyrmion phase, taking advantage of lessons gained in constructing suitable models in two-dimensions, which are discussed in Section II. We first explain why such a three-dimensional counterpart phase is anticipated, before explaining how toy models for the 3D phase are constructed taking inspiration from the topological phase known as the Hopf insulator 41 . To do this, we explain the defining characteristic of the 3D chiral skyrmion phase and how analogous physics defines the Hopf insulator phase, before utilizing techniques from Section I to construct a toy model Hamiltonian for the 3D chiral topological skyrmion phase from a two-band Hamiltonian for the Hopf insulator and a two-band Hamiltonian for its C ′ partner. We then introduce the topological invariant for characterization of the 3D chiral skyrmion phase, and discuss its numerical computation. We then use this invariant to characterize some canonical toy models for the 3D chiral skyrmion phase. Finally, motivated by past work on realization of unpaired Majorana zero-modes in two-dimensional systems with a defect that are characterized by the Hopf insulator Hamiltonian, we show that the Hamiltonian for the 3D chiral skyrmion phase with the same boundary conditions also realizes topologically-protected zero-modes.

B
Construction of models for the three-dimensional chiral skyrmion phase While a three-dimensional counterpart of the chiral topological skyrmion phase is expected in Hamiltonians with C ′ symmetry, as they correspond to the non-trivial homotopy group π 3 (Sp(2N )/U(N )) = Z 2 , it is not immediately clear how to construct a model realizing this topological phase. We expect the three-dimensional phase to correspond to formation of a hedgehog in the ground state spin expectation value over the three-dimensional Brillouin zone, rather than a baby skyrmion as in the two-dimensional case. Given that the atomic model for the two-dimensional chiral skyrmion phase is Eq. 7, which may be constructed from a two-band Chern insulator and its C ′ partner, we can see that the atomic model for the three-dimensional case would be constructed from some two-band Hamiltonian that realizes a momentumspace hedgehog and its C ′ partner. There is such a two-band model: Hamiltonians realizing the Hopf insulator 41 phase of matter realize such hedgehogs in the spin or pseudospin texture over the Brillouin zone when the topological invariant characterizing the phase, the Hopf invariant, is non-trivial. The two-band Bloch Hamiltonian h(k) for the Hopf insulator takes the following form: where for which we take z = (z 1 , z 2 ) T where The four-band atomic model for the three-dimensional topological skyrmion phase of matter may then be written as  Fig. 9 and 10, respectively. We show that the system realizes a skyrmion phase for a range of parameters ∆ 0 and M as characterized by the skyrmion number.

C Topological invariant in three-dimensions
Here, we introduce the topological invariant of the three dimensional topological skyrmion phases of matter. To provide context, we first review the Hopf invariant. Although there are some similarities between the form of the Hopf invariant and the topological invariant of the 3D topological skyrmion phase, they are computed with different physical quantities and are therefore in general distinct. The Hopf invariant is the topological charge of a vector field defined over the 3D Brillouin zone, where the vector is the projector onto occupied states. Following the paper introducing the Hopf insulator protected by C ′ symmetry in systems with 2 or more bands 45 , we consider a Hamiltonian H(k) which takes the following form when diagonalized, Continuously deforming all of the positive eigenvalues to +1, and all of the negative eigenvalues to −1, the deformed Hamiltonian takes the form, where U (k) is the matrix of eigenvectors and I N,N = I N ⊗ σ z , where I N is the N × N identity matrix and σ z is the Pauli matrix σ z = diag(1, −1). We can then define the Hopf invariant in terms of the projectors onto occupied states as in Liu et al. 45 as and therefore denote it as a projector topological invariant.
If the only degree of freedom labeling the fermionic creation and annihilation operators in a given k-sector is spin, this projector also corresponds to the ground-state spin expectation value. If the operators are labeled not just by spin, but also by other degrees of freedom, such as an orbital degree of freedom, or a particle-hole degree of freedom, the projector onto occupied states becomes distinct from the groundstate spin expectation value vector. In this case, one can still compute the topological charge for the projector vector field, but also a second topological charge, for the spin expectation value vector. In some cases, this spin topological charge quantizes and serves as a second topological invariant distinct from the projector-based topological invariant, which we define as a spin topological invariant. The 3D topological skyrmion phase corresponds to the system realizing non-trivial, quantized spin topological charge, when the spin expectation value vector is distinct from the projector onto occupied states. For these cases where the projector topological invariant and spin topological invariant are computed with distinct physical quantities, we may then write the topological invariant of the 3D chiral skyrmion phase as where here A Q is the skyrmion connection defined in Eq. 19, and Ω Q is the skyrmion curvature defined in Eq. 3. In this work, we restrict ourselves to four band models constructed from a two-band Hamiltonian and its C ′ partner, which each realize a topological phase characterized by the Hopf invariant. We may therefore compute a Hopf invariant using projectors onto occupied states as defined above, but also a separate 3D skyrmion number.
The 3D skyrmion number Q 3D may therefore be computed using Eq. 4 for the skyrmion curvature, and Eq. 20 for the skyrmion connection, computing the derivatives numerically. Using this invariant Q 3D , we characterize the topology of the toy models presented for the 3D chiral skyrmion phase.
A similar set of phase diagrams is computed for the Hamiltonian in Eq. (32) with free parameter M plus additional spin triplet pairing term Eq. (27) with free parameter ∆ 0 . The results are shown in Fig. 9 and 10. We show that the system can be in the three-dimensional topological skyrmion phase for a range of parameters ∆ 0 and M as characterized by the skyrmion number.

D Bulk-boundary correspondence
We now explore the bulk-boundary correspondence associated with non-trivial Q 3D . Bulk-boundary correspondence is expected when a projector topological invariant is non-trivial, and when symmetries protecting the topological phase in the bulk are respected in the bulk and by open boundary conditions 11,50 . Here, we instead investigate whether the bulkboundary correspondence is realized for a topological phase characterized by trivial projector topological invariant, the Hopf invariant n ω , and non-trivial spin topological invariant Q 3D .
To provide context, we first review the bulk-boundary correspondence of the Hopf insulator, as there are at least two different kinds of bulk-boundary correspondence, which have been considered for this system and are relevant to the present work. These two types of bulk-boundary correspondence are distinguished by different open boundary conditions. The first kind of bulk-boundary correspondence of the Hopf insulator is discussed in Moore et al. 41 introducing the Hopf insulator phase. This corresponds to open boundary conditions in thê z-direction, and periodic boundary conditions in thex andŷ directions. Under these boundary conditions, non-trivial Hopf invariant corresponds to gapless boundary modes known as Fermi rings localized on the top and bottom surfaces.
Here, it is important to distinguish between the Hopf insulator phase realized in a two-band system, and the Hopf insulator realized in systems with more than two bands protected by generalized particle-hole symmetry C ′43 . The twoband Hopf insulator has Z topological classification and Fermi ring surface states due to bulk-boundary correspondence. The C ′ -protected Hopf insulator with more than two bands instead has Z 2 topological classification. Furthermore, as the C ′ operation involves particle-hole conjugation in combina- 11. (a) The two-dimensional lattice with the defect located at (0, 0) and θ = arctan(y/x). (b) Schematic illustration of the annulus geometry considered for low-energy analytical theory of the defect bulk-boundary correspondence. For the Hopf insulator, at θ = 0, one unpaired Majorana zero-mode is localized on the inner boundary and one on the outer boundary of the annulus, respectively, as shown in this schematic. tion with spatial inversion, opening boundary conditions in theẑ-direction generically breaks C ′ near the surfaces. Gapless surface states are then not topologically-protected. As our 3D chiral skyrmion phase is protected by C ′ symmetry, we do not expect to find this kind of bulk-boundary correspondence.
A second kind of bulk-boundary correspondence has also been considered in the context of the Hopf insulator 51 . In this case, boundary conditions are first opened in thex-andŷdirections, while periodic boundary conditions are retained in theẑ-direction. The Hamiltonian is therefore defined in terms of two real-space coordinates x and y, and one momentum component k z . Then, one performs the substitution k z → θ(x, y). This spatially-varying θ is interpreted as an angle in the x − y plane, which can characterize a zero-dimensional defect.
We first review the results of Yan et al 44 on this second kind of bulk-boundary correspondence, in particular discussing the nature of the gapless boundary modes. The boundary conditions are shown in Fig. 11 (a) for the Hopf insulator lattice model. A corresponding annulus geometry for the lowenergy theory is shown in Fig. 11 (b). This geometry yields a bulk-boundary correspondence when the Hopf invariant is non-trivial, and the Hopf insulator Bloch Hamiltonian realizes a hedgehog in its (pseudo)spin texture centered at the origin of the three-dimensional Brillouin zone. According to the topological classification of defects by Teo and Kane 40 , the defect can correspond to a non-trivial homotopy group and realize a topologically non-trivial, odd-parity state. As the defect is a zero-dimensional system within the Brillouin zone, it is, however, not immediately clear how to realize unpaired Majorana zero-modes from this odd-parity state. A key insight of Yan et al. 44 was that re-interpreting the Hopf insulator Bloch Hamiltonian as describing a two-dimensional system with a defect means the zero-dimensional defect of the three-dimensional Brillouin zone corresponds to an effectively one-dimensional system along the line θ = 0 in the 2D system with a defect. Then, similarly to one-dimensional wire proposals for realizing unpaired Majorana zero-modes, an unpaired Majorana zero-mode is realized at each end of this θ = 0 wire. For such a finite-size, two-dimensional system described by the Hopf insulator Hamiltonian with k z replaced by angle θ characterizing a defect at the center of the system, we first numerically compute the spectrum and eigenstates of the Hopf insulator Bloch Hamiltonian as is done in Yan et al 44 , with a defect located at (0, 0) and θ = arctan(y/x). The results are shown in Fig. 12 for a finite-size system with 40 × 40 lattice sites. These results include the low-energy spectrum, the particle and hole component probability densities of a zeroenergy eigenstate, and the probability densities for each unpaired Majorana zero-mode, shown in Fig. 12 (a), (b), (c), and (d), respectively. As expected, the unpaired Majorana zero-mode at the boundary occurs for θ = 0, corresponding to the location of the zero-dimensional defect in the threedimensional Brillouin zone of the Hopf insulator.
We then additionally characterize the resultant unpaired Majorana zero-modes in greater detail. Specifically, after taking linear combinations of the two lowest energy eigenstates of the Hopf insulator Hamiltonian with these boundary conditions to compute arrays with 40 × 40 × 2 entries representing the unpaired Majorana zero-mode state wavefunction, we then plot each of these entries in the complex plane as a point, to better understand later results for the 3D chiral skyrmion phase. These representations in the complex plane of the unpaired Majorana bound states localized at the defect and boundary are shown in Fig. 12 (e) and (f), respectively. The MZM at the boundary appears as a straight line in the complex plane. This is expected for an unpaired MZM, which has a purely real wavefunction up to some phase degree of freedom. The structure of the MZM at the defect, in contrast, is polluted. This complicated structure of the zero-mode state localized at the defect is unexpected for an unpaired MZM: a key feature of unpaired MZMs is that they should be straight lines i.e. they are purely real up to a phase degree of freedom 52,53 . The unpaired MZM localized at the defect likely loses this feature because of coupling to higher energy states at θ ̸ = 0 due to the small radius of the defect.
Given that the 3D chiral skyrmion phase toy model presented here is constructed from the Hopf insulator Hamiltonian, we anticipate a similar bulk-boundary correspondence for the 3D chiral skyrmion phase, when describing a twodimensional system with open boundary conditions and a defect, as occurs for the Hopf insulator Hamiltonian itself. We consider the following 3D chiral skyrmion phase Hamiltonian, where H(k) is the 4-band Hopf Hamiltonian in Eq. (32) and ∆ t (k) is the spin triplet pairing matrix in Eq. (27). For this model, we then consider open boundary conditions in thex-andŷ-directions, and substitute an angle θ for momentum component k z . The results for 3D chiral skyrmion phase Hamiltonian Eq. (37) under these conditions, when the 3D skyrmion number is 1, are shown in Fig. 13 with additional disorder-averaging to illustrate that the results are independent of the specific form of ∆ t (k). The on-site disorder term takes the form of ϵ ∼ τ 3 σ 0 · ϵ 0 U(−1, 1) where τ 3 σ 0 is a C ′ conserving term, ϵ 0 is the magnitude of the disorder which is 5% of the difference between non-zero energy levels (i.e. the lowest positive level and the highest negative level), and U(−1, 1) is a uniform distribution with values from -1 to 1. The results shown are obtained from averaging 100 realizations.
While two zero-energy states are realized for the Hopf insulator under these boundary conditions as shown in Fig. 12(a), four zero-energy states are realized in the lowenergy spectrum under the same conditions for the 3D topological skyrmion phase, as shown in Fig. 13(a). The robustness of these states against local disorder in the latter case has been numerically verified by introducing random on-site perturbation at each lattice site. Fig. 13(b) shows the probability density profile of the particle component (denoted by u) and hole component (denoted by v) of a zero-energy eigenstate realized by the 3D topological skyrmion phase Hamiltonian, respectively. This can be compared directly with Fig. 12(b). Notably, we see that, unlike in the case of the Hopf insulator, the zero-modes of the 3D topological skyrmion phase outputted by numerical diagonalization are particle-hole symmetric, but also localized. Fig. 13(c)(d) each show one localized zero-mode for direct comparison with Fig. 12(c)(d). These two localized zeromode wavefunctions at the defect and boundary are then visualized in the complex plane in Fig. 13(e) and (f), respectively, similarly to Fig. 12(e) and (f). That is, each wavefunction is represented by an array of 40 × 40 × 4 complex numbers, and these complex numbers are plotted in the complex plane as individual points. Similarly to the unpaired MZM at the defect in the Hopf insulator, the zero-modes localized at the defect in the 3D topological skyrmion phase have complex structure. The zero-modes of the 3D topological skyrmion phase localized at the boundary, however, reveal an interesting difference between the case of the Hopf insulator and the 3D topological skyrmion phase. They appear as crosses, or two straight, perpendicular lines, rather than as a single straight line as for a MZM, in the complex plane. In agreement with the lowenergy theory presented in the following section and in the supplemental materials for the zero-modes on the boundary, linear combinations of these two-cross states yield two unpaired Majorana zero-modes (one purely real and one purely imaginary), each localized at the boundary. Rather than combine to form a generic complex fermion bound state, however, they are forbidden from hybridizing by a combination of C ′ and C symmetry to instead yield the cross structure.

E Analytical characterization of topologicallyprotected zero-modes
To further characterize the cross-zero modes of the 3D topological skyrmion phase, we derive a low-energy effective theory for the model Hamiltonian Eq. (37) with open boundary conditions corresponding to an annulus geometry as in Yan et al. 44 and shown schematically in Fig. 11. Full details are presented in the Supplementary Materials, and here we focus on the resultant low-energy effective theory. For a semi-infinite geometry, with the sample occupying the x > 0 region the resultant effective Hamiltonian written in the on-site basis of the full Hamiltonian (generalized particle-hole DOF and spin DOF) is By taking k y → −i∂ y and λ → y/R, H eff becomes Squaring both sides of the time-independent Schrödinger equation yields the resultant equation where ϵ y,± = −∂ 2 y + 1 R 2 y 2 ± 1 R + 1 4 ∆ 2 0 . The eigenstates of this low-energy effective Hamiltonian may then be expressed as {1, 2, 3, 4}) are eigenvectors of τ 2 σ 2 , and |ψ i ⟩ = |ψ i (y)⟩ is a spatially-varying scalar function satisfying the i th differential equation of the set of four differential equations contained in the Hamiltonian. This yields |µ i ⟩ which take the forms 41) These states serve as zero-energy eigenstates of the effective low-energy Hamiltonian if the following differential equations are satisfied by the {ψ i (y)} (omitting 1 4 ∆ 2 0 as higher order terms): Note that by setting ∆ 0 = 0 we recover the differential equation in the two-band Hopf insulator which only has a solution for − 1 R , so that we can take this solution, treat −i∂ y ∆ 0 as a perturbation term, apply first order perturbation theory (for degenerate energy levels), and obtain the spatial distribution of the zero-energy modes, |ψ 1 (y)⟩.
The zero-energy eigenstates at the inner boundary then take the form, The above analysis is carried out for x > 0 i.e. the inner boundary. When instead considering the other boundary corresponding to x < 0, we reach the same differential equations (S62) and (S64), this time with eigenstates |µ 2 ⟩ and |µ 4 ⟩ in Eqn (S61). Hence, we similarly find two additional zeroenergy states of the form, where |ψ 2 (y)⟩ is a corresponding spatial distribution for states at the outer boundary.
These results indicate that C ′ symmetry alone yields two Majorana zero-modes at the defect and two at the boundary, but this is insufficient to realize the cross zero-modes. C ′ symmetry in combination with C symmetry of the Hopf insulator Bloch Hamiltonian is also required: the cross states are the most general complex states permitted by C and C ′ symmetry. Specifically, certain linear combinations of these four cross states satisfy constraints imposed by C, that the four zeromodes may be considered as two pairs of MZMs, with one MZM in each pair localized at the defect and the other localized at the boundary, but may also instead be considered as two other pairings of MZMs as required by C ′ symmetry, with one pair localized at the boundary and one pair localized at the defect. This difference stems from the relationship between C and C ′ , that C ′ acts like C upon exchange of position and momentum variables. This point was, to our knowledge, first observed in work by Liu et al. 45 : under C transformation, the position variables are invariant, but all of the momenta variables pick up a sign. However, under the C ′ transformation, the momenta do not change, but all of the position variables pick up a sign. Therefore, for all intents and purposes, we can replace C ′ with C as long as we reverse the roles of the position and momenta variables.
We may also understand the non-trivial topological state from the perspective of the homotopy group π 3 (Sp(2N )/U(N )) = Z 2 being relevant both to mappings from the Brillouin zone to the space of flat-band Hamiltonians, but also to mappings from the Brillouin zone to the space of ground state spin expectation value yielding the topological momentum-space, ground-state spin textures discussed here. This could correspond to a Z 2 × Z 2 topological classification for the system overall. As a result, even if the topological invariant associated with mappings from the Brillouin zone to the space of flat-band Hamiltonians is even and trivial, the system may still harbor topologically-protected zero-modes due to the co-existing topological skyrmion phase.
This Hamiltonian realizing the cross zero-modes may therefore be interpreted as a generalized BdG Hamiltonian for a centrosymmetric superconductor, where C and C ′ do not transform the normal state Bloch Hamiltonian in the same way because it is not even in momentum. For such cases, the Hilbert space is enlarged with both a particle-hole DOF, and an additional generalized particle-hole DOF for C ′ . We can therefore also interpret the full Hamiltonian as a spinless superconductor, written in terms of the Bloch Hamiltonian as Here, c k is the annihilation operator of electrons with momentum k = (k x , k y , k z ).
It is interesting that the cross zero-modes are strongly localized only at one boundary or the other: this is due to only two linearly-independent cross states being possible at a given point in real-space. The cross-states therefore must additionally be orthogonal as a result of their spatial distribution. Therefore two cross states are strongly localized at the inner boundary, and the remaining two are strongly localized at the outer boundary to form a four-fold degenerate ground state manifold. Importantly, this indicates the cross zero-modes should be more robust against finite-size effects than unpaired Majorana zero-modes. This is supported by numerics, which show well-defined cross zero-modes at the outer boundary in smaller systems (see 13 (f)).
In summary, in the topologically non-trivial regime of the Hamiltonian presented here realizing the 3D topological skyrmion phase, corresponding to a skyrmion forming in the ground-state spin expectation value texture at half-filling, there are two cross states localized on the outer boundary, and two other cross states localized at the defect. These cross states are a manifestation of a bulk-boundary correspondence of topological skyrmion phases due to topology of defects trapped by the momentum-space spin textures.

IV Discussion & Conclusion
In this work, we present toy models for recently-identified topological phases of matter in two dimensions, characterized by the topological charge of a momentum-space spin texture, a spin topological invariant rather than the topological charge of a momentum-space texture of the projector onto occupied states, a projector topological invariant, as done in the ten-fold way classification scheme 54 . We realize these topological phases in systems which can therefore be more fully characterized in terms of both a projector topological invariant and a spin topological invariant. Extending our results to introduce three-dimensional topological skyrmion phases corresponding to trivial projector topological invariant and non-trivial spin topological invariant, we furthermore find evidence of a bulk-boundary correspondence of topological skyrmion phases. As bulk-boundary correspondence is generally expected only for non-trivial projector topological invariant, these results serve as an important signature and consequence of topological skyrmion phases.
The bulk-boundary correspondence considered here, a defect bulk-boundary correspondence, is realized by opening boundary conditions of the 3D system in two directions x and y, and replacing the momentum component in a third z direction by an angle θ(x, y), which can characterize a zerodimensional defect in the xy-plane. For these open boundary conditions, we characterize the system's topology in terms of a bulk projector invariant, the Hopf invariant, and a bulk spin invariant, the 3D skyrmion number. Under these open boundary conditions, the topological skyrmion phases corresponding to trivial projector invariant and non-trivial spin invariant possess topologically-protected "cross" zero-modes, which generalize unpaired Majorana zero-modes: the zeromode wavefunctions in lattice-based systems consist of a set of complex numbers that form a cross-like structure when each of these complex numbers is plotted in the complex plane, rather than forming a straight line corresponding to a purely real state up to some phase degree of freedom as in the case of an unpaired Majorana zero-mode.
In order to demonstrate this defect bulk-boundary correspondence of the topological skyrmion phases, we first present a method for constructing Bloch Hamiltonian models for topological skyrmion phases with four bands and characterize myriad examples. We furthermore show BdG Hamiltonians for centrosymmetric superconductors naturally possess the C ′ symmetry in addition to particle-hole symmetry C when the normal state Hamiltonian is even in momentum. More generally, the topological skyrmion phases are described by generalized BdG Hamiltonians corresponding to artificially enlarging the Hilbert space to possess a generalized particle-hole DOF in addition to a particle-hole DOF.
We introduce an efficient method for computing skyrmion numbers based on past work by Fukui et al. 46 . At the fourband level, we find the skyrmion number computed using the ground state spin expectation value at half-filling is no longer equal to the total Chern number. Instead, we find that the total Chern number C and the skyrmion number Q are related by C = −2Q analytically for the atomic model-and numerically for all other models considered-of the chiral topological skyrmion phase characterized by Bloch Hamiltonians with four bands, when each of C and Q are topologically stable. It is also possible for C to be topologically stable and Q to be topologically unstable, due to the magnitude of the ground state spin being zero somewhere in the Brillouin zone. We do not find any instances where the reverse case occurs, in which Q is topologically stable and C is topologically unstable, and will explore this further in future work.
We do not explore bulk-boundary correspondence of twodimensional topological skyrmion phases in this work. However, since completion of this work, bulk-boundary correspondence for the 2D topological skyrmion phases has also been further characterized in terms of the observable-enriched partial trace and observable-enriched entanglement spectrum introduced by Winter et al. 35,55 . That is, bulk-boundary correspondence associated with non-trivial skyrmion number Q, when boundary conditions are opened in one direction, corresponds to Q chiral modes in the observable-enriched entanglement spectrum bulk gap.
We use knowledge from constructing various toy models for the chiral topological skyrmion phase in two dimensions to introduce the three-dimensional topological skyrmion phase, constructing it from the Hopf insulator Hamiltonian. We characterize the topology of these models using both a projector topological invariant and a spin topological invariant. We then show that the three-dimensional spin skyrmion in momentum space traps a defect, which can yield a bulk-boundary correspondence for the three-dimensional skyrmion phase: taking a momentum component to be an angle θ characterizing a realspace defect in a two-dimensional system with open boundary conditions in the other two directions, we find two zeromodes are localized at the defect and two are localized at the boundary when the system is in the three-dimensional topo-logical skyrmion phase. These zero-modes are realized when the projector topological invariant, the Hopf invariant, is trivial and the spin topological invariant, the 3D skyrmion number, is non-trivial. Therefore, while bulk-boundary correspondence has previously been associated strictly with non-trivial projector topological invariants, we show bulk-boundary correspondence is also associated with non-trival spin topological invariant. Such a zero-mode exhibits a cross-like structure when each of the complex numbers characterizing the zeromode wavefunction is plotted in the complex plane. Such a cross zero-mode is protected by a combination of particle-hole and generalized particle-hole symmetries and is robust against disorder. They generalize unpaired Majorana zero-modes, and may be thought of as an unpaired Majorana zero-mode and its C ′ partner combined in a single, localized state, yet isolated from one another by a π/2 phase difference in the complex plane.
In conclusion, these results indicate that the topological skyrmion phases exhibit a bulk-boundary correspondence by trapping defects that can yield topologically-protected, gapless boundary states. This serves as confirmation of the topological skyrmion phases and their implications, including their type-II topological phase transition contradicting the flat-band limit assumption. As the flat-band limit assumption is used in the construction of the ten-fold way topological classification scheme of Schnyder et al. 54 and Ryu et al. 50 , the type-II topological phase transition of the topological skyrmion phases 35 indicates topological skyrmion phases are a new topological class outside the classification of Schnyder et al. 54 and Ryu et al. 50 , rather than additional topological characterization of the current topological classification of free fermions. In addition, these gapless boundary states due to defect bulk-boundary correspondence of the topological skyrmion phases may interact with the boundary states of other topological phases. Such interplay between myriad topological phases will be explored in future work, especially given the potential for the topological skyrmion phase and the Chern insulator phase to co-exist. Finally, the properties and consequences of the cross zeromodes of the three-dimensional topological skyrmion phase will also be explored further in future work, given the significance of unpaired Majorana zero-modes for topological quantum computation and the fundamental generalization of unpaired Majorana zero-modes presented by the cross zeromodes discovered here.
Supplemental material for "Defect bulk-boundary correspondence of topological skyrmion phases of matter" S1 Low-energy, analytical theory of the 3D topological skyrmion phase Here we derive a low-energy theory for the four-band, three-dimensional skyrmion phase with Hamiltonian defined in Eq. (34), for a semi-infinite geometry with the sample occupying the x > 0 region. Translational symmetry in thex-direction is therefore broken while momentum component k y remains as a good quantum number. Thus, a real-space coordinate is used in thê x-direction to label lattice sites. The wave functions and energy eigenvalues can be obtained by solving the following equation where ψ j is the j th entry of an eigenstate for this equation corresponding to the j th lattice site in thex-direction. Each ψ j has four components corresponding to the generalized particle-hole and spin degrees of freedom at each site in real-space, and the non-zero terms of the Hamiltonian matrix representation on the left-hand side of Eq. S1 take the following forms: the term on the diagonal is and terms off the diagonal are and where λ is a parameter substituted for theẑ-component of the momentum, k z , characterizing a defect. Here, λ = nθ(i, j) with θ being the relative angle between the defect and site(i, j), n is an integer taken to be 1 in this study and ∆ 0 is the strength of the spin-triplet pairing.
which is identical to the i = 1 case and is unnormalizable.
The above analysis is carried out for x > 0 i.e. the inner boundary. As for x < 0 (the outer boundary), the essential difference is that we swap H 1 with H † 1 and H 2 with H † 2 in Eqn (S1): where ϵ y,± = −∂ 2 y + 1 R 2 y 2 ± 1 R . Subsequently we will reach the same differential equations (S62) and (S64), again with eigenstates |µ 1 ⟩ and |µ 3 ⟩ in Eqn (S61). Hence, we similarly find two additional zero-energy states of the form, Thus, if the system has C ′ symmetry alone, there will be no localized Majorana zero modes which are equal-weight superposition of particle and hole components.
In the case of the skyrmion phase in the main text, |0⟩ R(I) is further constrained by |0⟩ I(R) , as the linear combinations |0⟩ R(I) ± |0⟩ I(R) must also respect this particle-hole symmetry. These constraints are severe enough to enforce a cross-like structure on each of the four localized eigenstates of the Hamiltonian at approximately zero energy, such that each consists of two Majorana zero-modes which are now isolated from one another by a π/2 phase in the complex plane, rather than being isolated from one another by separation over real-space.