Excitonic skin effect

TL;DR

This paper introduces the excitonic skin effect, a boundary localization phenomenon for particle-hole pairs caused by strong interactions and imaginary vector potentials, with implications for non-Hermitian physics.

Contribution

It presents a novel boundary localization effect for excitons in non-Hermitian systems, linking strong interactions, imaginary vector potentials, and boundary phenomena.

Findings

Boundary localization of particle-hole pairs due to imaginary vector potentials

Exponential growth of pairing effects in a non-Hermitian bosonic Kitaev model

Potential for experimental realization in atomic and electronic systems

Abstract

We show that strong interactions combined with band-dependent imaginary vector potentials give rise to boundary localization of particle-hole pairs, which we term the excitonic skin effect. In a bilayer system with layer-specific gain/loss and an in-plane magnetic field, excitons experience a net imaginary vector potential, resulting in directional amplification of particle-hole pairs. Including nearest-neighbor interactions leads to a non-Hermitian bosonic Kitaev model, where the pairing effects grow exponentially with the size of the system, revealing a unique form of critical skin effect in interacting systems. Our framework applies to both atomic and electronic platforms and is directly testable in current experiments. These results also provide a route to explore non-Hermitian analogs of tensor gauge fields.

Equations

$$\begin{array}{rl}& {\hat{H}}_s=\hat{K}+{\hat{H}}_d,\\ & \hat{K}=(\sum _{j,\sigma }-t{e}^{i{\varphi }_{\sigma }}{\hat{c}}_{j,\sigma }^{†}{\hat{c}}_{j+1,\sigma }+\Omega \sum _j{\hat{c}}_{j,\uparrow }^{†}{\hat{c}}_{j,\downarrow }+h.c.),\\ & {\hat{H}}_d=i\gamma \sum _j({\hat{c}}_{j,\uparrow }^{†}{\hat{c}}_{j,\uparrow }-{\hat{c}}_{j,\downarrow }^{†}{\hat{c}}_{j,\downarrow }),\end{array}$$

$$\hat{h}_{k}=\left(\begin{array}[]{cc}-2t\cos[(k-k_{0})a]+i\gamma&\Omega\\ \Omega&-2t\cos[(k+k_{0})a]-i\gamma\\ \end{array}\right),$$

$$E_{\pm }(k)=\alpha (k-i{A}_{\pm }{)}^2+O(k^3{a}^3)+C_{\pm },$$

$${\hat{K}}_e=\sum _{j,s=\pm }-(t_{R,s}{\hat{c}}_{j+1,s}^{†}{\hat{c}}_{j,s}+t_{L,s}{\hat{c}}_{j,s}^{†}{\hat{c}}_{j+1,s})+\Omega ({\hat{n}}_{+}-{\hat{n}}_{-}),$$

$$A_{\pm }=\frac{1}{2a}\ln (t_{L,\pm }/t_{R,\pm })=\frac{\mp \gamma \tan \varphi }{a\Omega }+O(\gamma /\Omega {)}^3,$$

$$\begin{array}{rl}& {\hat{V}}_{on}=V_0\sum _j{\hat{n}}_{j,\uparrow }{\hat{n}}_{j,\downarrow },\\ & {\hat{V}}_{nn}=V_1\sum _j{\hat{n}}_{j\pm 1,\uparrow }{\hat{n}}_{j,\downarrow },\\ & {\hat{U}}_{nn}=U_0\sum _j({\hat{n}}_{j\pm 1,\uparrow }{\hat{n}}_{j,\uparrow }+{\hat{n}}_{j\pm 1,\downarrow }{\hat{n}}_{j,\downarrow }),\end{array}$$

$${\hat{V}}_{on}=V_0\sum _j{\hat{n}}_{j,+}{\hat{n}}_{j,-}.$$

$${\hat{H}}_{ex}=\sum _j{J}_R{\hat{d}}_{j+1}^{†}{\hat{d}}_j+J_L{\hat{d}}_{j-1}^{†}{\hat{d}}_j,$$

$$A_{ex}=\frac{1}{2a}\ln (J_{L,\pm }/J_{R,\pm }).$$

$$|\psi (0)⟩={\hat{c}}_{j_0,+}^{†}{\hat{c}}_{j_0,-}|v⟩,$$

$${\hat{H}}_{ring}\sim$$

$$+e^{A_{wx}a}{\hat{c}}_{j,+}^{†}{\hat{c}}_{j,-}{\hat{c}}_{j+1,-}^{†}{\hat{c}}_{j+1,+}.$$

$$\begin{array}{rl}{\hat{U}}_{nn}=& U_0\sum _j\{({\hat{n}}_{j,+}+{\hat{n}}_{j,-})({\hat{n}}_{j+1,+}+{\hat{n}}_{j+1,-})+\\ & ({\hat{d}}_j^{†}+{\hat{d}}_j)({\hat{d}}_{j+1}^{†}+{\hat{d}}_{j+1})\},\\ {\hat{V}}_{nn}=& \frac{V_1}{2}\sum _j\{({\hat{n}}_{j,+}+{\hat{n}}_{j,-})({\hat{n}}_{j+1,+}+{\hat{n}}_{j+1,-})-\\ & ({\hat{d}}_j^{†}+{\hat{d}}_j)({\hat{d}}_{j+1}^{†}+{\hat{d}}_{j+1})\}.\end{array}$$

$${\tilde{H}}_{ex}=$$

$$+\Delta ({\hat{d}}_j^{†}{\hat{d}}_{j+1}^{†}+{\hat{d}}_{j+1}{\hat{d}}_j),$$

$$\hat{H}_{f}=\frac{1}{2}\hat{\Psi}_{f}^{\dagger}\left(\begin{array}[]{cc}\hat{H}_{11}&\hat{H}_{12}\\ \hat{H}_{12}^{T}&-\hat{H}_{11}^{T}\end{array}\right)\hat{\Psi}_{f},$$

$${\text{Trid}[a,b,c]}=\left(\begin{array}[]{cccc}a&b&&\\ c&\ddots&\ddots&\\ &\ddots&\ddots&b\\ &&c&a\end{array}\right).$$

$$\{\begin{array}{l}\mathcal{P}={\mathcal{P}}_y={\sum }_j{\hat{n}}_{j,+}-{\hat{n}}_{j,-}\\ \mathcal{Q}={\mathcal{Q}}_{xy}={\sum }_{j}j({\hat{n}}_{j,+}-{\hat{n}}_{j,-}),\end{array}$$

$${\hat{c}}_{j,s}^{†}\to$$

$${\hat{c}}_{j,s}\to$$

$${\hat{K}}_e=\sum _{j,s=\pm }-\sqrt{t_{R,s}{t}_{L,s}}({\hat{a}}_{j+1,s}^{†}{\hat{a}}_{j,s}+h.c.)+\Omega ({\hat{n}}_{+}-{\hat{n}}_{-}).$$

$$\{\begin{array}{l}P=⟨\psi |{\sum }_j({\hat{\eta }}_{+}{\hat{n}}_{j,+}-{\hat{\eta }}_{-}{\hat{n}}_{j,-})|\psi ⟩\\ Q=⟨\psi |{\sum }_{j}j({\hat{\eta }}_{+}{\hat{n}}_{j,+}-{\hat{\eta }}_{-}{\hat{n}}_{j,-})|\psi ⟩.\end{array}$$

$${\mathcal{N}}_{ex}=⟨\psi |\sum _j{\hat{d}}_j^{†}{\hat{d}}_j|\psi ⟩.$$

$$N_{ex}=⟨\psi |\sum _j{e}^{2A_{ex}j}{\hat{d}}_j^{†}{\hat{d}}_j|\psi ⟩,$$

Full text

The skin effect—a phenomenon where an extensive number of physical states accumulate at system boundaries—has emerged as a striking feature in both quantum and classical systems Yao and Wang (2018); Kunst et al. (2018); Martinez Alvarez et al. (2018); Lee and Thomale (2019); Borgnia et al. (2020); Zhang et al. (2020); Okuma et al. (2020); Okuma and Sato (2023); Lin et al. (2023). In non-Hermitian systems, this concentration is induced by non-reciprocal couplings that push eigenstates toward the boundaries Nelson and Vinokur (1993); Hatano and Nelson (1996, 1997, 1998). Such a non-Hermitian skin effect has enabled advances in quantum sensing, photonics and wave control Regensburger et al. (2012); Feng et al. (2014); Wiersig (2014); Longhi et al. (2015); Wiersig (2016); Liu et al. (2016); Chen et al. (2017); Weidemann et al. (2020). In parallel, it has been found that curved spaces naturally support skin effect Zhang et al. (2021). Because of the geometric distortion inherent in curved manifolds, wavefunctions tend to localize and amplify in certain regions. Skin effects also arise in classical systems. For example, random walks on tree graphs exhibit directional drift and boundary accumulation Hughes and Sahimi (1982); Cassi (1990); Monthus and Texier (1996).

Although these skin effects were initially studied in isolation, recent work has revealed a deep connection between them Lv et al. (2022). Central to this unification is the concept of an imaginary vector potential Nelson and Vinokur (1993); Hatano and Nelson (1996). In non-Hermitian systems, it emerges from non-reciprocal couplings; in curved spaces, it arises naturally from the Laplace-Beltrami operator da Costa (1981); and in classical random walks, it manifests as a drift term. Unlike the ordinary real vector potentials, which imprint real phases to the wavefunctions, imaginary vector potentials induce position-dependent distortions, leading to the accumulation of physical states at boundaries, i.e., the skin effect.

While the skin effect has been extensively studied in non-interacting systems, recent research has begun to explore its manifestations in correlated systems, where interactions enrich the non-Hermitian dynamics Li et al. (2020); Yang et al. (2021); Shen and Lee (2022); Lee (2021); Kawabata et al. (2022); Alsallom et al. (2022). Notably, the concept has been extended to multipolar excitations, such as dipoles and quadrupoles, via non-reciprocal correlated hopping Gliozzi et al. (2024). These developments pushes the study of skin effects into a new regime, where composite excitations exhibit boundary accumulation, opening possibilities for higher-order topological control and interaction-enabled non-Hermitian devices.

Here, we point out a new class of skin effect—the excitonic skin effect—arising from the interplay between interactions and band-dependent imaginary vector potentials. A prototypical system exhibiting this behavior is a bi-layer structure subject to an in-plane magnetic field and the layer-dependent loss or gain, which generates a band structure where the ground and excited bands carry imaginary vector potentials of equal amplitude but opposite direction. When strong onsite interactions bind a hole in the ground band with a particle in the excited band to form an exciton, this composite experiences a finite net imaginary vector potential, leading to directional amplification and boundary accumulation of excitonic states—a hallmark of the excitonic skin effect.

Including nearest neighbor interactions further enriches the physics. The system then supports pair creation and annihilation processes, which map onto a non-Hermitian bosonic Kitaev model in the hard-core boson limit. The imaginary vector potential amplifies the pairing terms, with their influence growing exponentially with the system size. This behavior can be viewed as an interacting counterpart to the critical skin effect. In addition to electronic systems, our framework applies to spin-orbit coupled atoms under dissipation, making it also relevant for current ultracold atom experiments. It is also worth mentioning that viewing the energy as a synthetic dimension allows the band-dependent imaginary vector potential to be interpreted as a higher-rank imaginary tensor gauge field—an analog of real tensor gauge fields Pretko (2017); Ma et al. (2018); Bulmash and Barkeshli (2018); Zhang et al. (2025). Our work thus opens up a new venue for studying the non-Hermitian analogs of tensor gauge fields in terms of the excitonic skin effect.

We consider a bilayer system as shown in Fig. 1a. Two layers in the $x-y$ plane is separated from each other by a distance $d$ in the $z$-direction. A magnetic field is applied in the $y$-direction. Since the dynamics in the $y$-direction is regular and the skin effect exists only in the $x$-direction, we focus on the $x-z$ plane, which can be described by a two-leg ladder in Fig. 1b. We first start from the single-particle Hamiltonian, which is written as

[TABLE]

where $\sigma=(\uparrow,\downarrow)$ is a pseudospin index that denotes the upper and lower layer, respectively. $\hat{c}^{\dagger}_{j,\sigma}$ ($\hat{c}_{j,\sigma}$) is the fermionic creation (annihilation) operator at site $j$ for spin-$\sigma$, and $\hat{n}_{j,\sigma}=\hat{c}^{\dagger}_{j,\sigma}\hat{c}_{j,\sigma}$. $\hat{K}$ is a two-leg Harper-Hofstadter model, $\phi_{\uparrow}=-\phi_{\downarrow}=\phi$ is the Peierls phase, $2\phi$ is the magnetic flux per plaquette, and $\Omega$ is the intra-lyaer tunneling strength. $\hat{H}_{d}$ denotes the dissipation in the system. The upper and lower layers have gain and loss rate, $\pm\gamma$, respectively. It is worth mentioning that the dissipative rates in these two layers do not have to be equal in magnitude and opposite in sign. We note that $i\sum_{\sigma}\gamma_{\sigma}\hat{c}^{\dagger}_{j,\sigma}\hat{c}_{j,\sigma}=i\bar{\gamma}\sum_{\sigma}\hat{c}^{\dagger}_{j,\sigma}\hat{c}_{j,\sigma}+i\gamma(\hat{c}^{\dagger}_{j,\uparrow}\hat{c}_{j,\uparrow}-\hat{c}^{\dagger}_{j,\downarrow}\hat{c}_{j,\downarrow})$, where $\bar{\gamma}=(\gamma_{\uparrow}+\gamma_{\downarrow})/2$ and $\gamma=(\gamma_{\uparrow}-\gamma_{\downarrow})/2$. Any finite $\gamma_{\uparrow}-\gamma_{\downarrow}$ thus leads to $\hat{H}_{d}$. The term proportional to $\bar{\gamma}$ only leads to a trivial overall decay of the total particle number and does not change our main results of the skin effect.

$\hat{H}_{s}$ in Eq.(1) is directly realizable in ultracold atoms, where $\hat{H}_{s}$ corresponds to the Hamiltonian of spin-orbit coupled fermions with dissipation, $\sigma=\uparrow,\downarrow$ represent two hyperfine spin states, and $\Omega$ denotes the Raman coupling strength. The phase carried by Raman lasers is imprinted onto atoms such that a two-leg Harper-Hofstadter ladder with a finite magnetic flux $2\phi$ is accessed Lin et al. (2011); Galitski and Spielman (2013). $\hat{H}_{d}$ with a finite $\gamma$ has recently been introduced in experiments Ren et al. (2022); Tao et al. (2025). In solids, $\hat{H}_{s}$ is a standard Hamiltonian describing two coupled layers with a parallel magnetic field when $\gamma=0$. It applies to both ordinary bilayer systems or Moiré superlattices. A finite $\gamma$ can be accessed by introducing losses to one layer through photon or electric field induced emission of electrons. Alternatively, layer-dependent optical driving leading to excitations to higher energy bands could also effectively result in the desired losses.

A Fourier transform of $\hat{H}_{s}$ provides the Hamiltonian in the momentum space,

[TABLE]

where $k$ is the crystal momentum, $a$ is the lattice spacing in the $x$-direction, and $k_{0}=-\phi/a$ denotes the shift of the crystal momentum when the magnetic field is implemented. $\hat{h}_{k}$ can be easily diagonalized and the two eigenenergy bands $\{|+\rangle,|-\rangle\}$ have eigenenergies $E_{\pm}(k)$. When $\gamma\neq 0$, $E_{\pm}(k)$ becomes complex. A profound feature of $E_{\pm}(k)$ is that ${\text{R}e}E_{\pm}(k)\sim k^{2}$, ${\text{I}m}E_{\pm}(k)\sim k$ in a broad parameter regime. For instance, when $|\Omega|\gg|\gamma|$, $|\Omega|\gg|t|$ and $|ka|\ll 1$ are satisfied, $E_{\pm}(k)$ has a simple analytical form,

[TABLE]

where $\alpha=a^{2}t\cos(k_{0}a)$, $A_{\pm}=\frac{\pm\gamma\tan(k_{0}a)}{a\Omega}$, and $C_{\pm}$ are momentum-independent constants. A notable feature of $E_{\pm}(k)$ is that $A_{\pm}$ arise in the energy as imaginary vector potentials. Unlike a real vector potential that shifts the band structure in the real momentum axis, here, an imaginary vector potential provides a linear imaginary part of the energy and the real part remains quadratic, as shown in Fig. 1c.

A previous study has considered the imaginary vector potential in a single band and the subsequent single-particle skin effects Zhou et al. (2022). Here, we focus on a unique property of Eq.(3) that the ground and the excited bands are equipped with imaginary vector potentials of the same amplitudes and opposite directions. The origin of such a band-dependent imaginary vector potential can be unfolded if we rewrite Eq.(1) in a different basis. We define $\hat{c}_{j,\pm}=(\hat{c}_{j,\uparrow}\pm\hat{c}_{j,\downarrow})/\sqrt{2}$. In the limit $|\Omega|\gg|\gamma|$ and $|\Omega|\gg|t|$, a standard perturbation approach gives rise to an effective Hamiltonian,

[TABLE]

where $t_{R,+}=t_{L,-}=t(\cos\phi+\gamma\sin\phi/\Omega)$, $t_{L,+}=t_{R,-}=t(\cos\phi-\gamma\sin\phi/\Omega)$, and $\hat{n}_{s}=\sum_{j}\hat{c}_{j,s}^{\dagger}\hat{c}_{j,s}$. The effective Hamiltonian corresponds to two decoupled Hatano-Nelson chains with non-reciprocal tunneling $t_{L,s}\neq t_{R,s}$, as shown in Fig. 1d. Such non-reciprocal tunnelings correspond to imaginary vector potentials $A_{\pm}$, which can be extracted from $t_{R,\pm}=\bar{t}e^{-A_{\pm}a}$ and $t_{L,\pm}=\bar{t}e^{A_{\pm}a}$. Here $\bar{t}=\sqrt{\cos^{2}\phi-\gamma^{2}\sin^{2}\phi/\Omega^{2}}t$, and

[TABLE]

consistent with the result in the continuum limit in Eq.(3). $A_{\pm}$ give rise to the well-celebrated non-Hermitian skin effect. The two chains have an energy difference $2\Omega$, corresponding to the ground and excited bands. Since the imaginary vector potentials of these two bands point towards opposite directions, a band-dependent skin effect arises, and eigenstates in the ground and excited bands are localized at opposite edges of the systems.

We now take into account the interaction $\hat{H}_{int}$. While atoms, in general, have short-range interactions, some magnetic species or the Rydberg states exhibit dipole-dipole interactions, which may lead to considerable nearest-neighbor interactions Griesmaier et al. (2005); Stuhler et al. (2005); Gallagher and Pillet (2008); Browaeys et al. (2016). As for electrons, the generic expression for the Coulomb interaction $\hat{H}_{col}=\hat{V}+\hat{U}$ is written as $\hat{V}=\tilde{V}_{0}\sum_{j,m}(m^{2}+d^{2}/a^{2})^{-1/2}\hat{n}_{j+m,\uparrow}\hat{n}_{j,\downarrow}$, and $\hat{U}=\tilde{U}_{0}\sum_{j,m\neq 0}|m|^{-1}(\hat{n}_{j+m,\uparrow}\hat{n}_{j,\uparrow}+\hat{n}_{j+m,\downarrow}\hat{n}_{j,\downarrow})$. $\hat{U}$ is the intra-layer Coulomb interaction that decays as $1/m$ as the separation between two sites $m$ increases, and $\hat{V}$ is the inter-layer Coulomb interaction. In Moiré systems, the tunable lattice spacing of the superlattice provides experimentalists with an extra tuning knob to control $a/d$, the ratio of the inter-layer and intra-layer interaction. As such, without loss of generality, we consider the interaction $\hat{H}_{int}=\hat{V}_{on}+\hat{V}_{nn}+\hat{U}_{nn}$,

[TABLE]

where $V_{0}=\tilde{V}_{0}{a}/{d}$, $V_{1}=\tilde{V}_{0}(1+d^{2}/a^{2})^{-1/2}$ and $U_{0}=\tilde{U}_{0}$.

Since the interaction decays as the separation between two sites increases, we first consider the dominant term, the onsite inter-layer interaction, which is rewritten as

[TABLE]

The form for $\hat{V}_{on}$ remains unchanged as we rewrite it in terms of $\hat{n}_{j,\pm}=\hat{c}^{\dagger}_{j,\pm}\hat{c}_{j,\pm}$. We consider a particle-hole excitation on a fully filled ground band. In the limit of $|V_{0}|\gg|t|$, the particle in the excited band and the hole left in the ground band prefer to stay at the same lattice site to avoid the extra energy cost for a particle-hole pair to stay at different lattice sites. In other words, an exciton is formed. This exciton can tunnel in the lattice via a second-order process, as shown in Fig. 2a. The effective Hamiltonian for a single exciton is written as,

[TABLE]

where $\hat{d}_{j}=\hat{c}_{j,-}^{\dagger}\hat{c}_{j,+}$ and $J_{R(L)}=\frac{-2t^{2}}{V_{0}}(\cos\phi\pm\gamma\sin\phi/\Omega)^{2}$. We see that the tunneling of an exciton is also non-reciprocal since $J_{R}\neq J_{L}$. Such nonreciprocal tunneling implies that an exciton is subject to a finite imaginary vector potential $A_{ex}$. From $J_{R(L)}=\bar{J}e^{\mp A_{ex}a}$, we obtain

[TABLE]

A straightforward calculation shows that $A_{ex}=A_{+}-A_{-}$, where $A_{\pm}$ is the band-dependent imaginary vector potential of single particles in Eq. (5). This can be understood from the fact that a hole in the ground band is subject to an imaginary vector potential of the opposite direction compared to a particle in the same band, i.e., $-A_{-}$. The net imaginary vector potential applied to an exciton is thus the difference between $A_{+}$ and $A_{-}$.

A finite imaginary vector potential $A_{ex}$ acting on an exciton is expected to induce the excitonic skin effect. To confirm this, we numerically simulate the dynamics governed by the full model $\hat{H}_{s}+\hat{V}_{on}$. We consider the initial state

[TABLE]

where $|v\rangle=\prod_{j}\hat{c}_{j,-}^{\dagger}|0\rangle$ is an insulating state with the ground band fully occupied. $|\psi(0)\rangle$ thus corresponds to an excitation of a single exciton that occupies the $j_{0}$-th site. We numerically computed the time evolution of this initial state, $|\psi(t)\rangle=e^{-i(\hat{H}_{s}+\hat{V}_{on})t}|\psi(0)\rangle$, on a $7\times 2$ ladder where it has 7 sites along the x-direction, and 2 sites along the y-direction representing the ground and excited band. The exciton is prepared in the middle of the ladder at $j_{0}=4$.

Using exact diagonalization, we obtain the time dependence of the exciton density $n^{ex}_{j}(t)\equiv\langle\psi(t)|\hat{d}_{j}^{\dagger}\hat{d}_{j}|\psi(t)\rangle$ in Fig. 2b. When $t=0$, $n^{ex}_{j}=\delta_{j,j_{0}}$ . As time goes by, the exciton propagates in the lattice. A notable feature is that this propagation is not symmetric and gets amplified toward the right, and therefore the distribution of the exciton density $n^{ex}_{j}(t)$ has larger probability on the right-hand side. This is the characteristic feature of the skin effect in the presence of an imaginary vector potential $A_{ex}<0$ acting on the exciton. Unlike the boundary accumulation of single-particles in the ordinary skin effect, here, particle-hole pairs concentrate at the edge, signifying a new type of skin effect for excitons. It is worth pointing out that, in spite of the directional amplification, conserved quantities exist once we appropriately take into account the finite curvatures underlying Hatano-Nelson chains (supplementary materials).

Whereas the above discussions have clearly demonstrated the excitonic skin effect, it is useful to consider such results from a different perspective. In the language of synthetic dimensions, the energy could be regarded as an extra dimension $w$ perpendicular to the real dimensions. Eq. (4) thus includes an imaginary rank-2 tensor gauge field, $A_{wx}=\partial_{w}A_{x}=A_{+}-A_{-}$. Recent studies of fracton phases of matter Yoshida (2013); Vijay et al. (2015, 2016); Pretko (2017); Ma et al. (2018); Prem et al. (2018); Bulmash and Barkeshli (2018); Nandkishore and Hermele (2019); Pretko et al. (2020); Yuan et al. (2020); Lake et al. (2022) have shown that, when single-particle excitations become immobile, the couplings between higher-rank gauge fields with dipoles and other multipoles become critical. However, those studies have mainly focused on real tensor gauge fields. Here, an imaginary tensor gauge field $A_{wx}$ provides excitons with a complex Peierls phase, leading to an effective non-Hermitian ring exchange interaction,

[TABLE]

The excitonic skin effect thus could also be regarded as a consequence of an imaginary rank-2 tensor gauge field.

Now we consider effects beyond the on-site interaction. To this end, we express $\hat{H}_{int}$ in terms of $\hat{c}_{j,\pm}$ and $\hat{c}_{j,\pm}^{\dagger}$, then the nearest-neighbor interactions can be written as

[TABLE]

Consider the constraint that the total number of fermions per site is fixed at 1, $\hat{n}_{j,+}+\hat{n}_{j,-}$ is equivalent to an identity operator for any $j$. The first terms in the expressions for $\hat{U}_{nn}$ and $\hat{V}_{nn}$ can then be dropped off. We thus only need to concentrate on the terms that depend on the dipole operators $\hat{d}_{j}$ and $\hat{d}_{j}^{\dagger}$. These terms contain hopping, pair creation and annihilation of excitons on nearest-neighbor sites. The effective Hamiltonian for excitons is written as

[TABLE]

where $\Delta=U_{0}-V_{1}/2$, $\mu=2\Omega$, and $\tilde{J}_{R(L)}=\frac{2t^{2}d}{V_{0}a}(\cos\phi\pm\frac{\gamma\sin\phi}{\Omega})^{2}-\Delta$. When $\tilde{J}_{L}=\tilde{J}_{R}$, the above equation has the same form as the bosonic Kitaev model McDonald et al. (2018). Here, $\tilde{J}_{L}\neq\tilde{J}_{R}$, and each exciton should be regarded as a hard-core boson. $\tilde{H}_{ex}$ thus corresponds to a non-Hermitian bosonic Kitaev model in the hard-core limit. Alternatively, by mapping the hard-core bosons to spin-1/2s, Eq. (13) could be regarded as a non-Hermitian generalization of the spin Kitaev model.

A unique feature distinguishing $\tilde{H}_{ex}$ from its Hermitian counterparts is that the non-reciprocal tunneling significantly amplifies the pairing term. To show this drastic effect transparently, we solve the energy spectrum of $\tilde{H}_{ex}$ by mapping it to a fermionic model $\hat{H}_{f}$, which has the same energy spectrum as $\hat{H}_{ex}$, via a Jordan-Wigner transformation $\hat{d}_{j}^{\dagger}=\exp(-i\pi\sum_{n<j}\hat{f}_{n}^{\dagger}\hat{f}_{n})\hat{f}_{j}^{\dagger}$. The last two terms in Eq.(13) become the $p$-wave pairing for fermions, and $\hat{H}_{f}$ is a non-Hermitian generalization of the 1D Kitaev model Kitaev (2001). The conditions for the existence of Majorana zero modes in such a model have been previously discussed Li et al. (2022); Yan et al. (2023); Ardonne and Kurasov (2025). Here, we focus on how the non-reciprocal coupling and the pairing term control the energy spectrum.

$\hat{H}_{f}$ in the Nambu space is written as

[TABLE]

where $\hat{H}_{11}=\text{Trid}[\mu,-J_{L},-J_{R}]$, $\hat{H}_{12}=\text{Trid}[0,\Delta,-\Delta]$, and $\hat{\Psi}_{f}=(\hat{f}_{1},\cdots,\hat{f}_{L},\hat{f}_{1}^{\dagger},\cdots,\hat{f}_{L}^{\dagger})^{T}$. ${\text{Trid}[a,b,c]}$ is a triadiagonal matrix,

[TABLE]

We note that the two diagonal blocks can be interpreted as two Hatano-Nelson chains with opposite non-reciprocity, one for the particle and the other for the hole in the $p$-wave superconductor. As such, the pairing term $\Delta$ can be viewed as the diagonal couplings between nearest-neighbor sites of these two chains, as shown in Fig. 3a. Since this is similar to the critical non-Hermitian skin effect, where onsite couplings exist between two Hatano-Nelson chains with opposite non-reciprocity Li et al. (2020), it is desirable to explore whether the critical skin effect may exist in this interacting system.

We numerically solve the energy spectrum of $\hat{H}_{f}$ with the open boundary condition (OBC). As shown in Fig. 3c, for a finite system with a fixed number of sites $L$, when $\Delta=0$, all eigenenergies are real. Increasing $\Delta$ to a critical value $\Delta_{c}$, some eigenenergies become complex. In particular, the larger $L$ is, the easier it is for a system with a fixed $\Delta$ to acquire complex eigenenergies, as shown in Fig. 3b. In other words, $\Delta_{c}$ decreases exponentially with increasing $L$. The scaling behavior of $\Delta_{c}$ is shown in Fig. 3d, indicating an exponential enhancement of the pairing $\Delta$ as $L$ increases. We thus conclude that the nearest-neighbor interaction and the non-reciprocal tunneling give rise to an interacting counterpart of the critical non-Hermitian skin effect.

We have shown that the excitonic skin effect arises from the interplay of magnetic fields, interactions, and dissipation. This phenomenon is potentially observable in both atomic and electronic systems. Exploring this new type of skin effect may expand the frontier in the study of higher-rank tensor gauge fields to a complex domain. It will also provide a new platform for physicists to study mon-Hermitian physics of composite particles

The authors thank Ian Spielman, Emmanuel Gutierrez and Yihang Zeng for helpful discussions. This work is supported by The U.S. Department of Energy, Office of Science through the Quantum Science Center (QSC), a National Quantum Information Science Research Center, and the Air Force Office of Scientific Research under award number FA9550-23-1-0491.