High spin axion insulator

Axion insulators possess a quantized axion field θ = π protected by combined lattice and time-reversal symmetry, holding great potential for device applications in layertronics and quantum computing. Here, we propose a high-spin axion insulator (HSAI) defined in large spin-s representation, which maintains the same inherent symmetry but possesses a notable axion field θ = (s + 1/2)2π. Such distinct axion field is confirmed independently by the direct calculation of the axion term using hybrid Wannier functions, layer-resolved Chern numbers, as well as the topological magneto-electric effect. We show that the guaranteed gapless quasi-particle excitation is absent at the boundary of the HSAI despite its integer surface Chern number, hinting an unusual quantum anomaly violating the conventional bulk-boundary correspondence. Furthermore, we ascertain that the axion field θ can be precisely tuned through an external magnetic field, enabling the manipulation of bonded transport properties. The HSAI proposed here can be experimentally verified in ultra-cold atoms by the quantized non-reciprocal conductance or topological magnetoelectric response. Our work enriches the understanding of axion insulators in condensed matter physics, paving the way for future device applications.

Axion insulators possess a quantized axion field θ = π protected by combined lattice and time-reversal symmetry, holding great potential for device applications in layertronics and quantum computing.Here, we propose a high-spin axion insulator (HSAI) defined in large spin-s representation, which maintains the same inherent symmetry but possesses a notable axion field θ = (s + 1/2) 2 π.Such distinct axion field is confirmed independently by the direct calculation of the axion term using hybrid Wannier functions, layer-resolved Chern numbers, as well as the topological magneto-electric effect.We show that the guaranteed gapless quasi-particle excitation is absent at the boundary of the HSAI despite its integer surface Chern number, hinting an unusual quantum anomaly violating the conventional bulk-boundary correspondence.Furthermore, we ascertain that the axion field θ can be precisely tuned through an external magnetic field, enabling the manipulation of bonded transport properties.The HSAI proposed here can be experimentally verified in ultra-cold atoms by the quantized non-reciprocal conductance or topological magnetoelectric response.Our work enriches the understanding of axion insulators in condensed matter physics, paving the way for future device applications.
Symmetry plays an essential role in understanding the behavior of condensed materials [1][2][3][4] .For example, in the presence of time-reversal symmetry, three dimensional insulator typically falls into two categories: one is trivial insulator while the other is topological insulator 5,6 .These divergent categories can be well described within the framework of the Chern-Simons theory, where the Lagrangian incorporates an additional symmetry allowed term L θ = dtdr 3 αθ/(4π 2 )E • B with E and B the conventional electric and magnetic fields, α the fine structure constant, and θ the gauge dependent axion term 7 .Because of the 2π periodicity under a gauge transformation 8 , the axion term here is well defined within the region [0, 2π).Besides, as the quantity E • B flips sign under time-reversal (T ) operation, the axion field manifests only two distinct values, that is, θ = 0 for normal insulator and θ = π for topological insulator 7 .Furthermore, the non-vanishing axion term in the Lagrangian would introduce additional magneto-electric responses to the Maxwell equations and in turn, results in a distinctive topological magneto-electric effect [9][10][11] , furnishing a hallmark to the quantized axion field.
In addition to the time-reversal (T ) symmetry, the quantized axion field θ = π can also be protected by combined lattice and time-reversal symmetry (for example S = T τ 1/2 with τ 1/2 the half translation operator) 12 , as the quantity E • B undergoes the same sign reversal.This kind of materials, termed axion insulator [13][14][15][16][17][18][19][20][21][22] , holds significant potential in layertronics [23][24][25][26][27] and quantum computing 28,29 .MnBi 2 Te 4 and its family have recently been proposed as axion insulators in the antiferromagnetic state 13,[30][31][32][33][34] , which finds a concise description with an effective Hamiltonian defined on the orbital and spin-1/2 spaces 31 .Because the symmetry transformations of high spin representations and spin-1/2 are identical (see Supplementary Section 1), in this work, we generalize this model to other spin species and thus propose a high spin axion insulator (HSAI) preserving the same symmetry.We find that the HSAI with spin-s possesses a notable axion field θ HSAI = (s + 1/2) 2 π.It carries a multiple quantized helical hinge current (QHHC) that is robust against disorders even in the absence of the topologically protected gapless exciations, which contradicts the integer surface Chern number.Consequently, HSAI exhibits an unusual quantum anomaly that violates the conventional bulk-boundary correspondence.In contrast to the case of spin-1/2 axion insulator, the direct calculation of the axion term shows that the large axion field in high-spin case originates mostly from localized surface Wannier functions while, in the bulk, the axion field is either 0 or π.Strikingly, we show that the axion field in HSAI can be tuned precisely by manipulating the magnetic configuration through an external magnetic field, providing a pioneering tuning knob to control the QHHC and the quantized magneto-electric response.Possible experimental realization in ultra-cold atoms is also discussed.

Results
Effective model for the HSAI Recalling the effective four-band Hamiltonian for the spin-1/2 axion insulator 31 , we consider a generic model defined on the high spin space which can be expressed as Here, d 0,1,2,3 = m 0 − Bk 2 , Ak x , Ak y , Ak z with A, B, m 0 the system parameters.Γ 0 = s 0 ⊗ τ 3 and Γ i=1,2,3 = s i ⊗ τ 1 where s i and τ i are matrices defined on the high spin space and 2 × 2 orbital space, respectively.The momentum k = k x , k y , k z is defined on a cubic lattice with the lattice constant a 0 inside the first Brillouin zone.This model Hamiltonian is given directly from the spin-1/2 axion insulator.A construction from symmetry perspective is provided in Supplementary Section 2. It is evident that the first term in Eq. ( 1) describes a high-spin topological insulator, which preserves both time-reversal (T ) and parity (P) symmetries.The second term describes the exchange interaction between topological electrons and normalized magnetic spins m s = m x s , m y s , m z s with ∆ the exchange strength, resembling that in MnBi 2 Te 4 , hence explicitly breaks the time-reversal symmetry while preserves the S symmetry in the infinite size limit along z-direction.We consider the antiferromagnetic phase of an even-layer slab involving only the antiparallel (or canted) spins in the top and bottom layers as illustrated in Fig. 1a, which restores the combined parity and time-reversal (PT ) symmetry.Unless otherwise specified, we adopt the typical model parameters as follows: Figure 1(b) displays the two dimensional energy spectra of the spin-3/2 HSAI in the absence (blue dashed lines) and presence (red solid lines) of the magnetic exchange term.In the former case, the time-reversal symmetry is present, where the energy spectrum is gapped in the bulk but has two conducting surface states on each side (blue dashed lines).These two surface states can be accurately fitted by a massless Dirac band E 1 ∼ k and a cubic band E 2 ∼ k 3 (inset in Fig. 1b).Turning on the exchange term in the latter case opens a band gap as indicated by the red solid lines in Fig. 1b.In both cases, the energy spectra are doubly degenerated because of the inherent (T or PT ) symmetry.
Layer-resolved Chern number, quantized helical hinge current and quantum anomaly To explore the topological properties of the HSAI, we calculate the layer-resolved Chern number C z along ẑ-direction 30,35 along with the cumulative Chern number Cz = z −Lz/2 C z .Given the Chern number C = (s + 1/2) 2 in the odd-layer case 36 , the system is a high Chern number insulator as shown in Supplementary Section 6.In the even layer system, the opposite layer-resolved Chern numbers are overall confined antisymmetrically inside few surface (top and bottom) layers as shown in Fig. 1c, resulting in a vanishing total Chern number C = 0. Nevertheless, the surface Chern number on one side turns out to be well quantized [C

top(bot) surf
= ∓2] when s = 3/2 as long as the Fermi level lies inside the energy gap (Fig. 1d).Because the layer-resolved Chern number is related to the axion field through the relation θ HSAI = (C bot surf − C top surf )π 37 .The oppositely quantized surface Chern numbers in spin-3/2 HSAI thereby assure a quadruple axion field θ HSAI = 4π.Moreover, since the Chern number difference between neighboring top (bottom) and side surfaces is an integer, HSAI supports a possible hinge state that is absent in spin-1/2 axion insulator 24 , which allows a subsequent QHHC owing to opposite chiralities on different surfaces.Nonetheless, we find that this QHHC survives counterintuitively without the existence of any hinge state.
To clarify this, let us first examine the average position ⟨z/L z ⟩ on the front surface of the slab at y/L y = −1/2.The results shown in Fig. 2b reveals two branches (red lines in Fig. 2b) concentrated unexpectedly around the bottom hinge (z/L z = −1/2) and propagating rightwards because of the positive group velocity.By contrast, we observe two other branches concentrate oppositely around the top hinge (z/L z = 1/2), which, simultaneously, propagate leftwards as depicted by the blue lines.Note that only the results for the front surface (at y/L y = −1/2) are presented here.In the presence of PT symmetry, the energy spectrum in Fig. 2b is doubly degenerated as stated above.There are four additional branches existing on the other surface at y/L y = 1/2.Because the wavefunctions on the diagonal hinges are connected by this PT symmetry, they must propagate along the same direction, supporting a helical hinge current.
We then turn to the spectrum density A(k x , E) on a semi-infinite slab [38][39][40] , where the system extends infinitely along +ẑ-direction but remains unchanged along the lateral directions.A(k x , E) on the front lower hinge illustrated by the blue dashed line in Fig. 2a are plotted in Fig. 2c.It shows that A(k x , E) originates mostly from the right-moving energy bands, agreeing remarkably well with the average position in Fig. 2b.This spectrum density can be verified experimentally by using the nano angle-resolved photoemission spectroscopy and microscopy [41][42][43] .Besides, the high spectrum density on the hinge indicates the presence of a hinge current, the current density of which can be quantitatively determined by 24,40 where E F is the Fermi energy labeled by the white line in Fig. 2c, r = (y, z), H HSAI (k x , r) is the Hamiltonian for the HSAI and G r (E F ; k x , r) is the retarded Green's function.
The upper panel in Fig. 2d presents the hinge current density J x (z) = 0 y=−Ly/2 j x (r) as a function of the layer index z.We see that J x (z) is verily confined on the hinge, in agreement with ⟨z/L z ⟩ and A(k x , E). Strikingly, this hinge current decays oscillatively into the side surface, exhibiting a beating mode (magenta line) in sharp contrast to that in spin-1/2 axion insulator 24 .This peculiar behavior can be quantitatively fitted by the superposition of two power-law decaying edge currents 44 , where k F 1 and k F 2 are the Fermi momenta for the two distinct modes marked by the white stars in Fig. 2c, while a 1(2) and ϕ 1(2) are fitting parameters.The integral of the current density over the layer index provides the current flux I x (z) = z 0 dzJ x (z) (middle panel in Fig. 2d), which oscillates around 2e/h and coincides perfectly with the fitting data.Additionally, the z-averaged current ⟨I x (z)⟩ = z 0 dzI x (z)/z (red line) quantizes to 2e/h only a few layers away from the hinge.Imposing a finite thickness along ẑ-direction enables us to calculate the moving average current The result displayed in the bottom panel demonstrates a helical hinge current quantized precisely to ±2e/h.Although the HSAI supports a QHHC identical to its integer surface Chern number, the topologically protected gapless exciations in lower dimension are completely absent as plotted in Fig. 2b.It worths note that the energy gap in Fig. 2b may be induced by the finite size effect.Our analysis in Supplementary Section 3 shows that the size dependence of the energy gap comes from the bulk bands, which demonstrates that this energy gap originates from the magnetic exchange interaction.Consequently, the conventional bulk-boundary correspondence that an integer Chern number must hold chiral edge states fails in HSAI, establishing an unusual quantum anomaly.

Non-reciprocal conductance
Owing to the chirality bonded to the quantized surface Chern number, this QHHC can be unveiled by the non-reciprocal conductance 2e, where G ij is the differential conductance 24,45 .In this device, two longitudinal leads (terminals 5 and 6) are intimately connected to the two ends of the sample while four transverse leads (terminals 1, 2, 3, and 4) are attached to different hinges on the front surface.Figure .2f shows three representative non-reciprocal conductances versus the Fermi energy E F in the clean limit (solid lines), with non-magnetic Anderson disorder (dashed lines) and with magnetic Anderson disorder (dashed dotted lines).In general, G N 65 = 0, G N 31 = −2e 2 /h and G N 42 = 2e 2 /h when the Fermi level lies inside the band gap for all three cases, consistent with the current distribution in Fig. 2d as well as the layer-resolved Chern numbers in Fig. 1d.This verifies that the QHHC is topological protected as the quantized non-reciprocal conductance is immune to both non-magnetic and magnetic Anderson disorders that even breaks the global PT symmetry 5,6 .
Since multiple frequency ac current is robust against ambient perturbation, to detect this QHHC, we employ an alternate detection in which terminals 2, 4, 5, and 6 are grounded whereas a harmonic voltage V (t) = V 0 sin (ω 0 t) with a periodicity T = 2π/ω 0 is applied alternatively to terminal 1 or 3 as illustrated in Fig. 2g.During the first (second) half period, a positive (negative) voltage is applied to terminal 1 (3) as an input while the current flows i 3 (t) [i 1 (t)] from terminal 3 (1) is detected as an output.Their combination gives rise to an asymmetric net current i(t) = i 1 (t)+i 3 (t) as shown in Fig. 2h.Performing a Fourier transform converts i(t) into I(ω).The non-reciprocal conductance can then be determined from the equation with N the number of periodicity (see Supplementary Section 4 for details).The result displayed in Fig. 2i shows that F (ω) = G N 13 when ω = 2ω 0 .Thus, non-reciprocal conductance measurements offer a reliable experimental method to visualize the QHHC in HSAI.

Axion term
The HSAI can alternatively be characterized by the axion term 7 , which can be evaluated directly from the hybrid Wannier functions (HWFs) constructed in terms of the Bloch wavefunctions 46 .In this scenario, the axion term on a slab is defined as 46 where z nk is the hybrid Wannier charge center along ẑ-direction and Ωxy knn is corresponding non-Abelian Berry curvature.
Figure 3(a) shows z nk in the first Brillouin zone for a six-layer HSAI with spin-3/2.These z nk consist of two different types, those localized on the top and bottom surfaces as emphasized by the red and blue lines and those extending into the bulk denoted by black lines.Those surface Wannier bands will disappear under a periodic boundary condition when connecting the top and bottom surfaces.The total axion term of the slab can subsequently be divided into two parts θ slab CS = θ bulk CS + θ surf CS with θ bulk CS and θ surf CS the axion terms corresponding to the surface and bulk HWFs.The bulk axion term θ bulk CS is identical to that obtained analytically from the Chern-Simons three form in the infinite size limit 46 .In Fig. 3b, we show θ bulk CS (red upside down triangle), θ surf CS (black circle) together with θ slab CS (blue triangle) versus the inverse thickness 1/L z .There are three distinctive features in this figure.First, the total axion term shows an obvious tendency quantized to θ slab CS = 4π when the system size approaches infinity (1/L z → 0), which confirms the quadruple axion term in two dimensional HSAI slab.Second, the axion term originates completely from the surface HWFs although the bulk HWFs also result in a small value that decreases θ bulk CS at finite size.Third, the axion term θ surf CS obeys the relation HW F is the top (bottom) surface Chern number obtained from the surface HWFs as indicated by cyan diamond (magenta square).These peculiar results reaffirm the unusual quantum anomaly and also the quadruple axion term θ HSAI = 4π in HSAI with spin-3/2.
Topological magneto-electric effect Such quadruple axion term implies a unique topological magneto-electric effect 13,30 .When applying a magnetic field B z to the HSAI along ẑ-direction, the Hamiltonian in Eq. 1 acquires a Peierls phase 30 , which redistributes the electron charge Q(z) in accordance to the confined layer-resolved Chern number C z as shown in Fig. 3c.The ensuing charge polarization P = Lz/2 z=−Lz/2 zQ(z)/L z almost quantizes to P ≈ 4αϕ, where α is the fine structure constant and ϕ is the total magnetic flux penetrating the HSAI slab.In comparison, a quantized orbital magnetization can emerge under an external electric field E z when incurring a potential drop δU = eE z L z in the HSAI Hamiltonian 47 .The red square in Fig. 3(d) shows the orbital magnetization M as a function of δU , which agrees quantitatively well with the ideal case benchmarked by the black line.These two results independently certify the quadruple axion term θ HSAI in spin-3/2 HSAI.The slight deviation from the exactly quadruple value originates from the finite size effect, which is further revealed by the size scalings of the axion term θ slab CS /π, polarization coefficient P/(αϕ), and magnetization coefficient M/(αδU ) against the inverse layer thickness 1/L z in Fig. 3e.We also evaluate the axion term and the surface Chern numbers for spin-5/2 HSAI in terms of the HWFs (Fig. 3f).The results displayed in Fig. 3g

Tunable topological phase transition
In the presence of an in-plane magnetic field, the antiparallel magnetic moments in the top and bottom layers become canted with the canting angle γ proportional to the magnetic field strength as illustrated in Fig. 4a.In this case, the quantized axion field in the infinite size limit is protected by m x P symmetry where m x is the mirror plane normal to x-direction.In Fig. 4b, we compare two dimensional band gaps as functions of γ for spin-1/2 axion insulator and spin-3/2 HSAI.Because the exchange gap is determined by the perpendicular magnetization M z , the band gap for spin-1/2 axion insulator decreases monotonically as γ is enlarged and finally becomes zero when γ = π/2.On the contrary, the band for spin-3/2 HSAI exhibits a gap close at γ = π/4 as shown in Fig. 4c, suggesting a possible topological phase transition.Indeed, Figure 4e shows that the surface Chern number obtained using both the Bloch wavefunctions and the HWFs(Fig.4d) changes from +2 (−2) to +1 (−1) when γ = π/4.At this point, the HWFs are connected at the Γ point (Fig. 4d), therefore the Berry curvature and the surface Chern number can transfer from one side to the other 48 , leading to an axionic phase transition.Such topological phase transition is further affirmed by the axion field, the polarization and magnetization coefficients shown in Fig. 4f.

Discussion
The device application of axion insulators requires the fine-control of the transport signals such as the magneto-electric response or the QHHC, which are identical to the axion field.In spin-1/2 axion insulators, the axion term cannot be tuned without disrupting the existing S symmetry or refabricating the setup 49 .Nevertheless, because different surface bands shown in Fig. 1b can be coupled via the in-plane exchange interaction (M x s x ⊗τ 0 ), an apparent topological phase transition between axion insulators with different axion fields occurs in the HSAI.Consequently, the axion term θ HSAI (in unit of π), hence the QHHC G N ij (in unit of e 2 /2h) and the magneto-electric effect P/(αB) [M/(αE)] in HSAI can be precisely adjusted from 4 to 2 via the application of an external magnetic field.Thus, our work opens up an exciting possibility for the groundbreaking advancement in the practical application of axion insulators.
In conclusion, we have proposed a HSAI defined on the high spin space and shown that this HSAI possesses a multiple axion field protected by the combined lattice and time-reversal symmetry.Notably, the axion term in the bulk of a HSAI still quantizes to θ = 0 or θ = π while the surface of HSAI possesses a large axion term and a consistent integer Chern number, which can be tuned by manipulating the magnetic configuration through an external magnetic field.These results extend the scope of recently discovered axion insulator in magnetic topological materials.In ultra-cold fermions on a honeycomb lattice, the exchange gap can be introduced by complex nextnearest-neighbour tunneling terms through circular modulation of the lattice position 50 .We thus propose that our theory can be tested in high spin ultra-cold fermions on a stacked honeycomb lattice, where the non-reciprocal conductance can be detected by the orthogonal drifts analogous to a Hall current under a constant force to the atoms 50,51 .

Methods
Caltulations of the layer-resolved Chern number, magnetization, and polarization.In a HSAI slab with periodic boundary conditions along the lateral dimensions, the momenta k x and k y are good quantum numbers because of the translation symmetry.Therefore, the layer-resolved Chern number can be calculated by projecting the total Chern number into specific layer, which can be written as Here, E F is the Fermi energy, E m(n) (k) is the eigenenergy of H HSAI with |m k ⟩ (|n k ⟩) the corresponding eigenstates, Pz = |ψ z ⟩⟨ψ z | is the projecting operator.The integral is performed inside the first Brillouin zone.
Under an electric field E z along ẑ-direction, a potential drop occurs inside the HSAI slab along the same direction.The onsite energy in each layer acquires an additional value eE z z with z the layer index and the total potential drop in the HSAI slab is δU = eE z L z .The orbital magnetization can then be obtained accordingly by using where HHSAI = H HSAI + eE z z with Ẽm(n) and | m⟩ (|ñ⟩) its eigenenergy and eigenstate, respectively.
Applying a magnetic field B z along ẑ-direction introduces a gauge potential to the HSAI lattice and thus breaks the in-plane translation symmetry.Inside each unit cell, HSAI acquires a gauge field ϕ 0 = drA • r/Ψ 0 with Ψ 0 = h/(2e) the magnetic flux quantum.The total magnetic flux penetrating the HSAI slab is ϕ = B z L x L y .We adopt the Landau gauge A = −yB z , 0, 0 , so the translation symmetry along ŷ-direction is broken while that along x-direction sustains.In this case, the charge distribution induced by the magnetic field can be obtained by using the Green's function method, yielding where r = x, y, z and the Green's function G r (E, r) = (E + iη − H HSAI ) −1 with η the imaginary line width function.On the other hand, as k x is still a good quantum number, the charge distribution along ẑ-direction can be alternatively obtained by using Moreover, because only the negative charge originating from electrons are considered here in Eqs.6 and 7, to derive the unbalanced charge distribution and in turn the polarization, the uniform background charge compensating the positive ions in the lattice has to be removed from the results, which has the form q back = − Lz/2 z=−Lz/2 q(z)/L z because of the charge conservation.As a result, the charge distribution has the form Q(z) = q(z) + q back .The charge polarization can finally be expressed as Caltulations of the axion term using the hybrid Wannier function.In a HSAI slab, the hybrid Wannier wavefunction |h n,k ⟩ can be constructed from the Bloch wave function.We thus have |h n,k ⟩ = 1/2π π −π dk z |n k ⟩e −i(k•r+kzz) .In this case, the hybrid Wannier charge center takes the form z n k = ⟨h n,k |z|h n,k ⟩ 46 .To calculate the non-Abelian Berry curvature, we divide the twodimensional Brillouin zone into a regular mesh with b x and b y being the primitive reciprocal vectors that define the mesh.Then the gauge covariant Berry curvature has the form 52

Ωxy
where The wavefunctions constructed by a linear combination of the occupied bands at neighboring mesh point are | hn, Green's function method for calculating the differential conductance G ij .The differential conductance G ij corresponds to the transmission coefficient T ij from terminal j to terminal i, which can be derived by using the non-equilibrium Green's function method.Based on the Landauer-Büttiker formula 45 , the transmission coefficient T ij can be expressed as Fermi energy, η the imaginary line width function and Σ i the self energy due to the coupling to terminal i.To incorporate the disorders, we generate random potentials δE ∈ −W/2, W/2 for the non-magnetic Anderson disorders or δM z ∈ −W z /2, W z /2 for magnetic Anderson disorders at each site r, then add these random potentials to the Hamiltonian in the Green's functions.The results in the presence of disorders are calculated under 10 times average (Fig. 2f).

Spin matrices for different spins
The spin matrices for different spin species can be generated mathematically from the Clifford algebra [1, 2], which are spins

Hamiltonian for high spin axion insulator from symmetry perspective
The model Hamiltonian for the high spin axion insulator can be built from the high spin topological insulator preserving both parity and time-reversal symmetry [3].These parity and time-reversal symmetry on the orbital and spin basis are defined as P = τ z and T = e −isy K, respectively, where τ z is the Pauli matrix acting on orbitals, s y is the spin matrix shown in Sec. 1, and K is the complex conjugation operator.Since (the odd power of) the spin matrices s x,y,z flip sign under the time-reversal symmetry, the minimum Hamiltonian for a high-spin topological insulator preserving both parity and time-reversal symmetry thus reads where k 2 = k 2 x + k 2 y + k 2 z , m 0 , B, A i=x,y,z are system parameters.The first term in Supplementary Eq. ( 1) is the kinetic energy while the second term represents the spin-orbital coupling.Due to the time-reversal symmetry, the axion field of H 0 is quantized.In the main text, we take the spin orbit coupling A x = A y = A z for convenience.However, our theory is universal and is not limited to special parameters.It also worths note that the spin orbital coupling is crucial for the quantized axion field, which otherwise becomes θ 0 = 0 if A i=x,y,z = 0.In the presence of magnetic ordering, the time-reversal symmetry is explicitly broken.However, in the antiferromagnetic phase, the combined lattice and time-reversal symmetry is still well preserved because the magnetic moments on the adjacent layers are antiparallel.The magnetic layers introduces an exchange interaction into the system, which recasts the Hamiltonian as where ∆ is the exchange gap between the magnetic moment m s and the topological electrons with spin s.In three dimension limit, the system preserves the symmetry S = T τ 1/2 with τ 1/2 the half translation operator along z-direction.Whereas, this S symmetry breaks into combined parity and time-reversal symmetry that can be defined as PT = σ z e −isy K on a slab geometry, where σ z is the Pauli matrix switching the magnetic moments between the top and bottom layers

Finite size effect
To rule out the possibility of the finite-size effect induced energy gap in HSAI, we study the size dependence of the energy gap.Supplementary figures 1 a and b show the band structures of a one-dimensional nanowire with size L y = 20 and L z = 10 for the system in the antiferromagnetic and ferromagnetic cases, respectively.In the ferromagnetic case, the system is a spin-3/2 Chern insulator with a Chern number C = 4 (see Sec. 6).Therefore, there are four gapless edge states denoted by the red lines in Supplementary Fig. 1 b, in sharp contrast to the gapped band structure for a HSAI in Supplementary Fig. 1 a.On the other hand, figure 1 c shows the band gap for the HSAI as a function of L y .It is apparent that this gap is independent of the size L y .To further rule out the possibility of the finite-size effect along z-direction, we compare the band gap for the HSAI with the bulk band gap of the Chern insulator after removing the edge bands (red lines in Supplementary Fig. 1b).
The results plotted in Supplementary Fig. 1d show that the two band gaps as functions of L z coincide quantitatively with each other, which demonstrates that the band gap in HSAI is also induced by the bulk bands rather than the finite size effect due to the overlapping between edge states on top and bottom surfaces.

Alternate detection
The temporal dependent asymmetric current output in the main text has the form for an arbitary integer n, where G ij is the conductance from terminal j to terminal i as explained in the main text, V 0 and ω 0 are the amplitude and frequency of the alternative harmonic voltage input, respectively.Performing a Fourier transform with respect to t converts the current into the frequency domain, which yields Therefore, the non-reciprocal conductance can be unveiled by the function (solid magenta lines in the top and middle panels) coincide remarkably well with the fitting data (blue dots) obtained from the superposition of three power law decaying edge currents where k F i refers to the Fermi momenta marked by the white dots in Supplementary Fig. 3c, a i and ϕ 0 are fitting parameters, certifing that the beating mode of J x (z) comes from the superposition of the coherent edge currents on the same hinge.This beating mode occurs only in the HSAI since it originates from the superposition of two or more coherent edge currents while spin-1/2 axion insulator harbors only one edge current.On the other hand, one can see that ⟨I x (z)⟩ and ⟨I x (z)⟩ M A quantize to ±4.5e/honly a few layer away from the hinge.Thus, those signatures evidently confirms the quantized helical hinge current identical to the surface Chern number.and robust against both non-magnetic and magnetic disorders.As a result, the quantized helical hinge current is topological protected.This quantized helical hinge current can be there exist two regions (magenta regions) where the band gap vanishes.In theses regions the system is a metal, therefore the axion term is not well-defined [5].It is also important to note that the topological phase transition originates solely from the surface axion term θ surf CS while the bulk axion term θ bulk CS remains unchanged during the process.Because the transport signals in HSAI such as the non-reciprocal conductance identical to the quantized helical hinge currents and the topological magneto-electric response are proportional to the axion field, the HSAI thus provides a platform to realized the goal of axionic topological phase transition.
6 Even-odd effect of the high spin axion insulator In the odd layer system, the net magnetization is non-vanishing because of the uncompensated magnetic layer.We can thus simulate this state by using parallel magnetic moments on both top and bottom layers.In order to reveal the difference between them, we explore the layer resolved Chern numbers, which can be obtained by using Eq. ( 4) in the methods.
The results displayed in Supplementary Fig. 7 show that the layer-resolved Chern numbers for odd layer systems distribute symmetrically on the top and bottom layers, leading to a vanishing axion field θ = 0 while a nonzero Chern number C = (s + 1/2) 2 analogous to the

Figure 1 :
Figure 1: Model of the HSAI.a Schematic for the HSAI defined on the |s, m z ⟩ space.The blue arrows represent the electron spin with different magnetic quantum number m z , which takes values ranging from −s to s individually.b Energy spectra of the spin-3/2 HSAI along M→ Γ →R path on a slab of thickness L z with (red solid lines) and without (blue dashed lines) the magnetic exchange interaction.Here, the black lines refer to bulk bands.Inset: Energy dispersion for the spin-3/2 HSAI in the absence of magnetic exchange term near the charge neutral point (solid lines) and the fitting data (markers).c Layer-resolved Chern number C z and the cumulative Chern number Cz = z −Lz/2 C z versus the layer index z for the spin-3/2 HSAI.The surface Chern numbers C top(bot) surf that summarize the layer-resolved Chern number on the upper (lower) half side is −2 (+2).d Surface Chern number as a function of the Fermi energy E F for the spin-3/2 HSAI.Here, the thickness of the HSAI slab is L z = 20.

Figure 2 :Figure 3 :
Figure 2: Transport properties of the spin-3/2 HSAI. a Schematic current flow in a HSAI.The red arrows denote the quantized helical hinge currents.b Energy spectrum and the average position ⟨z/L z ⟩ on the front surface for a HSAI nanowire with L y = 30, L z = 16.c Spectrum density A(k x , E) for the front lower hinge as labeled by the blue line in a on the k x − E plane.Here, the system size is L y = 30, L z = ∞/2.The white dashed line represents the Fermi energy E F = 0.1.The white stars that mark the intersects between the Fermi energy and the spectrum are the Fermi momenta k F 1 and k F 2 .d Top and middle panels are the current density J x (z), current flux I x (z) and its z-averaged flux ⟨I x (z)⟩ versus the layer index z for a semi-infinite system with size L y = 30, L z = ∞/2.The blue dots are the fitting data.Bottom panel shows the distribution of the moving averaged current ⟨I x (z)⟩ M A on the front surface with system size L y = 30, L z = 150.e Bird's eye view (top panel) and high-angle shot (bottom panel) for the six terminal device.f Ensemble-averaged non-reciprocal conductances versus the Fermi energy in the clean limit (W=0), with non-magnetic Anderson disorders of strength W = 1 and with magnetic Anderson disorders of strength W z = 0.3.Here, the system size is L x = 31, L y = 20, L z = 21, and the size of transverse terminals is 10 × 10. g Experimental setup to detect the nonreciprocal conductance.In this setup, terminals 2, 4, 5, and 6 are grounded.The voltage is applied alternatively to terminal 1 or terminal 3. h Corresponding temporal dependent current output with parameters G 13 = 4.5e 2 /h, G 31 = 2.5e 2 /h.i F (ω) as a function of the frequency ω.

. 1 :
Size dependence of the energy gap. a and b show the band structures for systems in the antiferromagnetic and ferromagnetic cases, respectively.In the ferromagnetic case, the system is a spin-3/2 Chern insulator.The red lines in b denote the four gapless edge states.Here, the system size is L y = 20, L z = 10.c, energy gap for the HSAI versus L y .d, size dependence of the HSAI band gap (red triangle) and Chern insulator bulk band gap (blue up-side-down triangle) with L y = 20.The results for the Chern insulator bulk band gap is obtained after removing the surface states (red lines).All parameters for the model Hamiltonian are exactly the same as those in the main text.
with N the truncation of the summation n.If ω = 2ω 0 , F (ω) = G 13 − G 31 = G N 13 with G N 13 the non-reciprocal conductance when N approaches infinity.However, the function F (ω) = 0 if otherwise.Consequently, the quantized helical hinge current can be Supplementary Fig. 3: Quantized helical hinge current in HSAI with spin-5/2.a Schematic current flow in a spin-5/2 HSAI slab.b One dimensional energy spectrum and the average position ⟨z/L z ⟩ on the front surface as functions of k x with L y = 20, L z = 24.c Spectrum density A(k x , E) for the front lower hinge as marked by the purple star in a on the k x − E plane.Here, the system size is L y = 20, L z = ∞/2.d Top and middle panels are the current density J x (z), current flux I x (z) and its z-averaged flux ⟨I x (z)⟩ versus the layer index z for a semi-infinite system with size L y = 20, L z = ∞/2.The blue dots are the fitting data using three power law decaying edge currents with momenta k F 1 , k F 2 and k F 3 as marked by the white dots in c.Bottom panel shows the distribution of the moving averaged current ⟨I x (z)⟩ M A on the front surface with system size L y = 20, L z = 150 on the neighboring hinges.To further quantitatively uncover those helical hinge currents, we calculate the current density J x (z), current flux I x (z), z-average current ⟨I x (z)⟩ and the moving averaging current flux ⟨I x (z)⟩ M A in accordance with the formalism provided in the Methods.The results are plotted in Supplementary Fig. 3d.On one hand, J x (z) and I x (z)

. 4 :
Transport properties of the spin-5/2 HSAI in a six terminal device.a Ensemble-averaged non-reciprocal conductances versus the Fermi energy in the clean limit, with non-magnetic Anderson disorders W = 1 and with magnetic Anderson disorders W z = 0.3.Here, the system size is L x = 31, L y = 20, L z = 21, and the size of transverse terminals is 10 × 10. b Corresponding temporal dependent current output between terminals 1 and 3 with parameters G 13 = 2.5e 2 /h and G 31 = 7e 2 /h.c F (ω) as a function of the frequency ω.In the same six terminal device as schematically shown in Fig. 2e in the main text, we calculate the non-reciprocal conductance for the spin-5/2 HSAI by using the non-equilibrium Green's function method.Three representative non-reciprocal conductances G N 31 , G N 42 and G N 65 as functions of the Fermi energy E F are plotted in Supplementary Fig. 4a under different conditions (in the clean limit W = 0, with non-magnetic Anderson disorders W = 1.0 and with magnetic Anderson disorders W z = 0.3).We see that G N 31 = 4.5e 2 /h, G N 42 = −4.5e 2 /h and G N 65 = 0 when the Fermi energy lie inside the band gap regardless of the presence of disorders, confirming that those quantized non-reciprocal conductances are chiral

Supplementary Fig. 7 :
Even-odd effect of HSAI. a and b are layer-resolved Chern number for spin-3/2 and spin-5/2 HSAI with different layer thicknesses L z = 19 (red triangle) and L z = 20 (blue square).The Chern numbers and axion fields are labeled in the figure.
demonstrate that spin-5/2 HSAI possesses a surface Chern number Systematic results for the spin-5/2 HSAI are provided in Supplementary Section 5.The topological properties for HSAI with different spin species are summarized in Table.1, giving a distinct axion field θ = (s + 1/2) 2 π and C

Table 1 :
Axion terms and surface Chern numbers for axion insulators with different spins.