Chiral symmetry breaking and topological charge of zigzag graphene nanoribbons

Although the self-consistent Hartree-Fock fields break chiral symmetry, our work demonstrates that zigzag graphene nanoribbons maintain their status as short-range entangled symmetry-protected topological insulators. The relevant symmetry involves combined mirror and time-reversal operations. In un-doped ribbons displaying edge ferromagnetism, the band gap edge states with a topological charge form on the zigzag edges. An analysis of the anomalous continuity equation elucidates that this topological charge is induced by the gap term. In low-doped zigzag ribbons, where the ground state exhibits edge spin density waves, this topological charge appears as a nearly zero-energy edge mode.


Introduction
Non-interacting zigzag graphene nanoribbons [1,2] exhibit chiral symmetry and are topological insulators [3].However, this symmetry breaks down in the presence of electron interactions, leading to uncertainty regarding the definition of a topological charge.Our paper aims to elucidate the formation of a topological charge in disorder-free interacting zigzag graphene nanoribbons.We will demonstrate that in the presence of such interactions, the relevant symmetry combines mirror and time-reversal operations.The breakdown of chiral symmetry leads us to utilize axial current instead of ordinary current for an accurate definition of a topological charge, mirroring the approach by Goldstone and Wilczek [4] in the context of solitons in polyacetylene [5].This issue holds experimental significance due to the recent fabrication of atomically precise nanoribbons (Cai et al. [6], Kolmer et al. [7], Houtsma et al. [8]).
One way to define a topological charge is through the generation of mass (or gap), which breaks chiral symmetry [4].This concept can be elucidated using a onedimensional Dirac equation.In one spatial dimension (and generally for all odd spatial dimensions), the gapless (or massless in the particle physics context) Dirac fermions possess the property of chirality.The Lagrangian density of the massless Dirac fermion (see the Appendix for details and conventions) is invariant under two kinds of global gauge transformations: ψ → e iΛ ψ ordinary, ψ → e iγ 5 Λ ψ axial or chiral, (1) where ψ is a Dirac spinor, and γ 5 = γ 0 γ 1 .Chirality is defined to be the eigenvalue of γ 5 , which can take only ±1 values.According to classical Noether's theorem, these invariances imply the conservation of the corresponding currents (Fermi velocity may be multiplied to the spatial component): Owing to a special property of gamma matrices in two dimensions, the above two currents can be related by: where ϵ µν is the two-dimensional Levi-Civita symbol with ϵ 01 = 1.Expressing the Dirac spinor ψ in γ 5 -diagonal basis as (ψ R , ψ L ) t ('t' is transpose), the explicit expressions of currents are given by (note that j 1 = −j 1 , etc.): A mass term for Dirac fermion couples ψ R and ψ L .Two independent mass terms (in Hamiltonian) are allowed in the following way: The above mass terms break the chiral symmetry explicitly.From the viewpoint of a one-dimensional polyacetylene system, these masses correspond to the even and the odd components of dimerization.In the presence of the mass terms, the continuity equation for the axial current is modified as follows [4,9]: Eq.( 6) remains valid even if the mass parameters depend on coordinates.For slow adiabatic temporal changes, the time derivative of Eq.( 6) can be neglected, and we have (ρ = j 0 denotes charge density).∂ρ ∂x = M(x).
An abrupt step-like change of the mass terms can give rise to a solitonic [5,10] topological charge [4] accumulated at the soliton location, where the spectral weight of a soliton is equally divided between conduction and valence bands.More precisely, M(x) in Eq.( 6) corresponds to the topological charge [11,12] density after subtracting out the trivial vacuum contribution [4].Then, the integration of Eq.( 7) over some interval with endpoints (say, [−∞, ∞]) gives the topological charge Q top associated with the solitonic configuration of mass parameters The determination of this topological charge relies on the topological aspects of the mass parameters, defined as M(x) = ∂θ ∂x , where the order parameter θ(x) represents a kink: In polyacetylene θ(x) represents dimerization order [5].However, this topological charge is not necessarily quantized in general circumstances [4].Then, by Eq.( 8), the induced fermion number can be identified with the topological charge Q top of mass parameters, which is in turn determined by their boundary behavior [10].Our goal in this paper is to try to understand the formation of a topological charge in disorder-free interacting zigzag graphene nanoribbons in view of the above statement.First, let us give a brief introduction to non-interacting graphene nanoribbons.The band structure of a non-interacting zigzag graphene nanoribbon is shown in Fig. 1.Zigzag graphene nanoribbons [2,[13][14][15] consist of A and B sublattices, as seen in the top figure of Fig. 2. The A/B sublattice degree of freedom is referred to as the chirality.The chiral transformation [16,17] on electron operators is implemented by where a i and b i are the electron destruction operators of A and B sublattice sites which are similar to the Dirac fermion ψ R , ψ L under the action of γ 5 (see the Appendix).Under the transformation Eq.( 10), the nearest-neighbor hopping term between A-B sublattices of the tight-binding Hamiltonian of graphene nanoribbon without on-site interaction changes sign, and this is called the chiral symmetry.This chiral symmetry implies a particle-hole symmetry within the energy spectrum and ensures the existence of zeroenergy soliton [18,19] edge states.Thus, non-interacting graphene zigzag nanoribbons are deemed a symmetry-protected topological insulator [20,21].Now, the inclusion of on-site repulsion interactions among electrons explicitly breaks chiral symmetry, as the density operator does not alter its sign under the transformation described by Eq. (10).This interaction introduces an energy gap, ∆.Despite this energy gap, the persistence of edge states is evident, as depicted in Fig. 2. Furthermore, interacting zigzag graphene nanoribbons manifest as short-range entangled symmetry-protected topological insulators [21].(Here, we exclude consideration of disordered zigzag ribbons that exhibit long-range entangled [22,23] topological order [21].)Consequently, it is natural to inquire about the topological origin of these edge states in the presence of on-site interactions.(In this paper, we do not consider armchair ribbons since the physical origin of a gap may differ.)

↑↓ ↑↓
< l a t e x i t s h a 1 _ b a s e 6 4 = " m B I P u n x 0 w X m a 8 I D / S w Y K J y 5 4 N Z J P W T s n t e P r s 9 L V W u p q H l 4 Q A O 4 R h c u I A K 3 E A V a k D h A Z 7 h F d 6 s F + v d + r A + J 6 0 5 a z q z D z O w v n 4 B J E q f I w = = < / l a t e x i t > S < l a t e x i t s h a 1 _ b a s e 6 4 = " B r K l h h 0 R 3 4 r l I u s E Y 3 U q i E 4 t B l e 4 c 1 6 s d 6 t D + t z 0 p q z p j P 7 M A P r 6 x c G V p 8 R < / l a t e x i t > A < l a t e x i t s h a 1 _ b a s e 6 4 = " b j U R c s y d u V z x 9 L t n l e s F 5 l R f The band structure of a non-interacting zigzag nanoribbon displays nearly zero-energy states localized on the zigzag edges.Second : Solitonic nearly zeroenergy states appear near k = ±π/a 0 , comprising symmetric ψ B and antisymmetric ψ A states well-localized on the left and right zigzag edges (a 0 represents the unit cell length of the zigzag ribbon).
We aim to address the aforementioned questions by investigating the interacting nanoribbons within the framework of the Hartree-Fock (HF) approach and employing an anomalous continuity equation similar to Eq. ( 6).The self-consistent HF fields break chiral symmetry.However, the self-consistent fields of the opposite zigzag edges exhibit antisymmetry under mirror symmetry operation.The relevant symmetry combines mirror and time-reversal operations, enabling the formation of a topological charge induced by the mass on ferromagnetic zigzag edges in undoped ribbons.This mechanism is analogous to the equation in Eq.( 8) and is of a non-perturbative nature.This topological charge persists even in the low-doped region where the ground state displays edge spin density waves.It appears as a nearly zero energy edge mode.Our work demonstrates that interacting graphene nanoribbons retain their status as symmetryprotected topological insulators.
We believe that the conventional classification scheme of topological insulators [24,25] does not directly apply to zigzag graphene nanoribbons because interacting zigzag ribbons display gapful edge excitations, whereas the conventional classification applies to gapless excitations.Moreover, the quasi-onedimensionality of zigzag ribbons is important as it spatially separates the opposite zigzag edges of different chiralities.

Tight-binding Hubbard model of Nanoribbon
Let us consider the narrowest armchair nanoribbon with zigzag edges and the narrowest zigzag nanoribbon, as illustrated in the top figure of Fig. 2. The Hamiltonian of the tight-binding Hubbard model of these nanoribbons [26][27][28][29] consists of the hopping term with hopping amplitude t and the Hubbard term for on-site Coulomb repulsion U B l e 4 c 1 6 s d 6 t D + t z 0 p q z p j P 7 M A P r 6 x c G V p 8 R < / l a t e x i t > , where i denotes a site, and σ represents spin up or down.A edge sites are predominantly occupied by spin-down, while the B edge sites are predominantly occupied by the opposite spin.Therefore, in each sublattice, the total occupation numbers of spin-up and spin-down electrons are different.Second: Band structures of two degenerate ground states of a zigzag nanoribbon in the presence of on-site Coulomb repulsion.Zigzag edge states are found near the band gap edges, occurring at energies approximately E ≈ ±∆/2, with momenta around k ≈ ±π/a 0 .L and R stand for left and right edge localized states.Third: Zigzag edge states ψ L and ψ R , exclusively localized on either the left or right edge, result from the combination of symmetric ψ S and antisymmetric ψ A solitonic states, as depicted in Fig. 1.
The second quantized hopping Hamiltonian is (h.c. is the Hermitian conjugate, and σ =↑, ↓ is spin) where c R,σ is the electron destruction operator at site R with spin σ and < R, R ′ > denotes a nearest-neighbor pair.One of the pair belongs to A-sublattice, while the other is in B-sublattice.The Hubbard term is given by where n Rσ = c † Rσ c Rσ is the electron number operator at site R with spin σ.As discussed in Section 1, H t posssses the chiral symmetry, while the H U explicitly breaks it.
Applying HF mean-field theory to Eq.( 13) at half-filling, a quadratic mean-field Hamiltonian H HF can be obtained: where the gap parameter ∆ Rσ , which corresponds to the mass parameter of Eq.( 6), is defined by A constant term has been added to H HF to make the Fermi energy to be zero.The site gap ∆ Rσ associated with spin σ can be understood as the staggered potential energy for an electron with the opposite spin σ.This potential is staggered since it changes sign between A-and B-sites, as depicted in Fig. 4. It is important to note that the mean-field HF Hamiltonian can be split into the spin-up and the spin-down components [30]: The band structure without the Hubbard interaction H U taken into account was shown in Fig. 1.The nearly zero-energy edge states exist close to the Brillouin zone boundary, and they exhibit non-magnetic behavior (the net spin polarization of a zigzag edge is zero).The existence of these zero energy chiral edge states crucially rests on the chiral symmetry.However, a self-consistent solution of the HF mean-field Hamiltonian reveals that for ribbon widths shorter than 100 Å, the ground state is antiferromagnetic in the undoped case.Local density approximation calculations [27,31,32] also yield similar results.In this state, site spins of sublattice A (B) are spin-up (-down), and the state is doubly degenerate with an energy gap.The chiral edge states do exist, and these chiral edge states are now oppositely polarized (edge antiferromagnetism) [26].Each zigzag edge displays ferromagnetism but the opposite zigzag edges with different chirality are antiferromagnetically coupled (see the top figure of Fig. 2).The band structures are displayed in the second figure in Fig. 2. The HF results qualitatively agree with those of density matrix renormalization methods [29], both in doped and undoped zigzag ribbons.In the low-doped regime, it is found that edge spin density waves develop in contrast to the ferromagnetic edges observed in the undoped case.Thus, doping is a singular effect.A staggered potential may be applied to graphene nanoribbons [33].In this case, mirror and time reversal symmetries are broken, as is their combined operation.The HF method shows that the charge of each zigzag edge may assume non-integer values instead of integer values.

Symmetry of undoped rectangular ribbons
Ribbons exhibit a symmetry based on the combined operation of parity transformation (mirror reflection) and time reversal transformations, as described by the relationship: This relation is illustrated in Fig. 3.Note that R is the mirror symmetric site of R and σ is the opposite spin of σ.Here, the relevant mirror axis of the zigzag nanoribbon is the horizontal x-axis and that of the armchair nanoribbon the vertical y-axis through the center of the ribbon.This symmetry aligns with the existence of a spin-up edge state on one zigzag edge and the corresponding spin-down edge state on the opposite zigzag edge.
Then the definition of the gap parameter Eq.( 15) immediately implies that Note that Eq.( 18) is valid even for the low-doped case as long as Eq.( 17) is valid.It is important to note that this mirror symmetry is accompanied by spin reversal and that the symmetry in the occupation numbers involves a spin reversal but not in the site gap values.At half-filling (undoped case), the condition ⟨ n Rσ ⟩ + ⟨ n Rσ ⟩ = 1 being combined with Eq.( 17) yields which is demonstrated in Fig. 4.

Symmetry of low-doped zigzag ribbons
Let us consider zigzag ribbons in the low doping limit.It turns out that doping has a singular effect [29]: the grounds state displays an edge spin density wave, as shown in We will call these axes the mirror axes.Note that only up-spin gap parameters are indicated.
Fig. 5.This is in contrast to the uniform edge ferromagnetism of the half-filled case.Note that in the low-doped region, the electron spins are collinear to a very good approximation, as confirmed by the density matrix renormalization group approach in the matrix product representation.This implies that our Hamiltonian at half-filling may still be applicable [29].
To understand the role of the doping we write the filling condition as where δ is the doping level.Then we have Eq. (21) shows that the antisymmetry of the gap parameter ∆ Rσ observed at half-filling in Fig. 4 is disrupted when moving away from the half-filling state relative to the mirror axis.

Anomalous continuity equation
The total current operator of the tight-binding Hamiltonian Eq.( 11) is given by which is defined such that the ordinary continuity equation is satisfied.Remember, in one-dimensional Dirac fermions, the axial current can be defined, where the roles of of the upper edge is out of phase with that of the lower edge.The site R is parametrized by (k, l).charge density ρ and spatial charge current density j interchange between ordinary and axial current (as shown in Eq.( 2)).Additionally, the Frölich equation governing the (unpinned) sliding charge density wave [34] is expressed as: Here, κ represents compressibility, v F signifies Fermi velocity, and E denotes an electric field.In one-dimensional systems, Eq.( 23) essentially embodies the well-known axial anomaly [9].These insights lead us to consider the time derivative of the spatial current, as outlined in Eq.( 22).This approach is motivated by the pronounced manifestation of chiral (axial) symmetry breaking effects, notably evident in the continuity equation for the axial current (Eq.( 6)).
The time derivative can be found from the commutator with H HF Eq.( 14): This expression comprises two contributions: one from the HF decomposed Hubbard interaction and the other from the hopping term.Taking the HF ground state expectation value of Eq.( 24) and considering the temporal adiabatic limit, (−iℏ) ∂⟨Jtot⟩ ∂t ≈ 0 we find : Then, we explore the implications of Eq.( 25) in light of the (integrated) anomalous continuity equation (Eq.( 6)) and the accumulation of the topological charge (Eq.( 8)).
The computation of the commutators are straightforward.For each spin projection, we find and It's important to clarify that in Eq.( 27), R 1 and R ′ are not a nearest neighbor pair, see Fig. 6.To elaborate further, if we consider a fixed R ′ , then R must be a nearest neighbor in the other sublattice, thereby determining the bond vector R − R ′ .Subsequently, R 1 should be a nearest neighbor of R, placing it in the same sublattice as R ′ .Now, considering the proximity of , effectively representing two times the number operator at R ′ .However, it's crucial to note that R 1 might not necessarily be the same as R ′ ; in such cases, R 1 becomes a next-nearest neighbor site from R ′ while still being in the same sublattice.Thus, Eq.( 27) can be expressed as Some examples of the next-nearest bond terms can be visualized using Fig. 7, such as ⟨ c † l+1,3,↑ c l,1,↑ ⟩.We note that in the framework of the anomalous continuity equation treated here, no bonds longer than the next-nearest ones appear.
It is important to emphasize that these results describe two-dimensional electron motion, covering the charge movement both along and perpendicular to the ribbon.In the following section we apply the above results to the zigzag nanoribbons depicted in the top figure of Fig. 2.

Topological charge of zigzag edges
In consideration of the anomalous continuity equation of the axial current (Eq.( 6)), we must compute the expectation value of the gap term (Eq.( 26)) and the terms originating from the hopping Hamiltonian (Eq.( 28)) concerning the HF ground state.Recalling that the HF Hamiltonian is split into the spin-up and the spin-down components (see Eq.( 16)) the expectation values can be computed for each spin component separately.Now it is crucial to note that there is a vertical mirror axis through the center of the zigzag ribbon (see Fig. 4).Then the expectation values of the bond operators of Eqs.(26,28) should respect this mirror symmetry.We hasten to comment that this mirror symmetry operation does not involve spin reversal.However, the mirror symmetry employed in the discussion of the gap parameter in Section 3 does involve spin reversal.Therefore, these two mirror symmetry operations are distinct from each other.
This vertical mirror symmetry implies that the horizontal x-component of the expectation values of the operators of Eqs.(26,28) vanish identically owing to the oddness of the horizontal component of the bond vector R − R ′ of Eqs.(26,28) under the mirror reflection, so that the topological charges can be accumulated only along the vertical direction: the upper and the lower chiral zigzag edges.Below l, m denote the unit cell and the position within the unit cell, respectively (see Fig. 7).
For simplicity, we'll evaluate the axial current for the smallest-width zigzag nanoribbons.Subsequently, we will validate this result for larger widths.This model of zigzag nanoribbons will clearly reveal the main physics.Firstly, regarding the expectation value of the hopping Hamiltonian in Eq.( 28), it comprises two distinct components.Initially, we focus on the first term, the contribution stemming from the density operator terms.We observe that the contributions from m = 2 and m = 3 sites (refer to Fig. 7) vanish.This occurs due to the vertical components of the nearest bond vector emanating from the site, summing up to zero.Thus, only the upper and lower edge sites contribute to the summation for the vertical y-component in the first term of So far we have investigated the topological charges of zigzag ribbons.It turns out that an interacting rectangular armchair ribbon with two zigzag edges and two armchair edges may also support edge states on its two zigzag edges.As the ribbon width increases more topological edge states are formed in the gap, as shown in Fig. 9.When a staggered potential is present [33], mirror and time reversal symmetries are broken, as is their combined operation.The HF method shows that the charge of each zigzag edge may assume non-integer values, i.e., they are not quantized.Plot of zigzag edge state probability density of undoped rectangular ribbons.As the ribbon width increases more nearly energy states are created in the gap ∆/2 = 0.25t.Here U = 2t.
We investigate whether ribbons in the low doping limit, characterized by the edge spin density ground state [29], possess topological edge states.Our numerical studies reveal that within the low-doped region, the ribbon can host nearly zero-energy gap states with a topological charge, as illustrated in Fig. 10.However, it's important to note that in the undoped case, the edge states have energies around E ≈ ±∆/2, as depicted in Fig. 2.

Discussion and Summary
Our findings suggest that chiral symmetry breaks due to on-site repulsion among electrons.However, an abrupt spatial change in the gap at the zigzag edges enables the emergence of topological charges with energies around ±∆/2, leading to the formation of a symmetry-protected topological insulator.This symmetry is dependent on combined operations involving parity (mirror reflection) and time reversal transformations.The importance of this combined symmetry may be further illustrated by considering a staggered potential [33].In the presence of such a potential, both mirror and time reversal symmetries are broken, as is their combined operation.The HF method shows that the charge of each zigzag edge may assume non-quantized charges as opposed to integer charges.Are topological edge states robust in the presence of disorder?The introduction of disorder induces a topological phase transition from a symmetry-protected insulator to a topologically ordered insulator (see Refs. [23,30]).Disorder in zigzag ribbons is a singular perturbation that gives rise to instantons.It also leads to the emergence of new edge states exhibiting e/2 fractional charges, which are protected by topological order, despite the presence of broken symmetries.But these new edge states are mixed chiral edge states, i.e., they are bonding and antibonding states of the chiral edge states.So, in this sense, the topological edge states are robust in the presence of disorder.
The boundary states (edge states) of our system are gapful.Therefore, it is not clear to what extent the classification scheme of Refs.[24,25] can be applied to our system.Moreover, the quasi-one-dimensionality of zigzag ribbons is important as it spatially separates the opposite zigzag edges of different chiralities.Disorder couples these edges and produce e/2 fractional charges.We believe that zigzag graphene nanoribbons are outside the conventional classification scheme.It would be interesting to expand the classification scheme to zigzag ribbons.
In low-doped zigzag ribbons, the ground state forms an edge spin density wave rather than a ferromagnetic edge.This state supports nearly zero-energy zigzag edge states, distinct from the ±∆/2 energies, exhibiting a gap-induced topological charge.
We hope that our work will stimulate the study of atomically precise nanoribbons recently fabricated [6,7] from the perspective of symmetry-protected topological insulators.It would be intriguing to observe topological charges of nearly zero energy states in low-doped ribbons that exhibit the edge spin-density wave ground state.Scanning tunneling microscopy [35] could be employed for this purpose.
t e x i t s h a 1 _ b a s e 6 4 = " m B I P u n x 0 P g C T K C F p d y H g 4 / p 6

Figure 2 .
Figure 2. First: Two interacting rectangular nanoribbons are shown.A graphene zigzag nanoribbon (a) and an armchair nanoribbon (b) consist of two sublattices, A (red sites) and B (blue sites).In an zigzag (armchair) ribbon, the left and right carbon sites are part of the armchair (zigzag) edges.The net site spins of the zigzag edges are displayed.The spins are 'down' on the red 'A' sites and 'up' on the blue 'B' sites.Zigzag edges are antiferromagnetically coupled, while each of them is ferromagnetically ordered.Site spin values are determined by the occupation numbers n iσ :s iz = 1 2 [n i↑ − n i↓ ], where i denotes a site, and σ represents spin up or down.A edge sites are predominantly occupied by spin-down, while the B edge sites are predominantly occupied by the opposite spin.Therefore, in each sublattice, the total occupation numbers of spin-up and spin-down electrons are different.Second: Band structures of two degenerate ground states of a zigzag nanoribbon in the presence of on-site Coulomb repulsion.Zigzag edge states are found near the band gap edges, occurring at energies approximately E ≈ ±∆/2, with momenta around k ≈ ±π/a 0 .L and R stand for left and right edge localized states.Third: Zigzag edge states ψ L and ψ R , exclusively localized on either the left or right edge, result from the combination of symmetric ψ S and antisymmetric ψ A solitonic states, as depicted in Fig.1.

Figure 3 .
Figure 3.The values of n i,↑ (represented by vertical solid lines) and n i,↓ (shown as vertical dashed lines) within a unit cell of the shortest-width zigzag ribbon are depicted.Note that for each site ⟨n R,↑ ⟩ + ⟨n R,↓ ⟩ = 1.

Figure 5 .
Figure 5. Spin-up and spin-down profiles of ⟨n R,σ ⟩ at the zigzag edges of a doped zigzag ribbon with an additional δN = 12 electrons added to half-filling.Note that ⟨n R,↑ ⟩ + ⟨n R,↓ ⟩ ̸ = constatnt, i.e., the varies from site to site.The site spins sR,z = 1 2 [⟨n R,↑ ⟩ − ⟨n R,↓ ⟩]of the upper edge is out of phase with that of the lower edge.The site R is parametrized by (k, l).

< l a t e x i t s h a 1 _Figure 6 .
Figure 6.The possible locations of sites represented by R 1 are indicated by dashed lines.

Figure 7 .
Figure 7. Unit cells of a zigzag ribbon are labeled by index l.
Figure 9.Plot of zigzag edge state probability density of undoped rectangular ribbons.As the ribbon width increases more nearly energy states are created in the gap ∆/2 = 0.25t.Here U = 2t.

EFigure 10 .
Figure 10.Occupied spin-up zigzag edge states of a doped ribbon δN = 12.There are also similar spin-down states.Here U = 2t.
Magnitude of the self-consistent gap ∆ R↑ at each site.Mirror symmetry about the horizontal (vertical) axis is present only in the armchair (zigzag) case.