Real-space decomposition of $p$-wave Kitaev chain

The exact solution of a model Hamiltonian, especially for many-body system, plays an important role in physics and sometimes may open the door to the exploration of new frontiers in physics. One of common and efficient tools is the Fourier transformation, which decomposes a Hamiltonian into many commutative sub-Hamiltonians due to the translational symmetry. In contrast to the $k$-space decomposition, there exist real-space decompositions, such as the block-diagonalization based on the reflection symmetry and so on. In this work, we propose another decomposition for the many-body Hamiltonian, which bases on an extended Bogoliubov transformations in real space for spinless fermions. By taking this transformation, we find that the Hamiltonians of a class of Kitaev chains Kitaev of length $2N$ with zero chemical potential can be mapped to two commutative ones of Kitaev chains of length $N$. All the eigen states can be constructed by the ones of two sub-systems. It provides a clear physical picture and a simple way to solve the dynamic problem in a complicated system from the result for relative simple systems. We demonstrate the implications of this decomposition by a SSH Kitaev model, which consists of dimerized hopping and pairing terms and supports the rich quantum phases. The features of the system, including the groundstate topology and nonequilibrium dynamics from a trivial initial state, can be revealed directly from that of well-known simple Kitaev chains. Based on this connection, two types of BCS-pair order parameters are introduced to characterize the phase diagram, showing the ingredient of two different BCS pairing modes.

For instance, the topological index of the original model can be shown to be simply the sum of that of two sub-models. The ground state of an SSH Kitaev model has rich quantum phases with winding numbers $\mathcal{N=}0$, $1$, and $2$, respectively, while a simple sub-system has quantum phases with winding numbers $\mathcal{N=}0$ and $1$. The physical picture becomes clear in the framework of decomposition. Based on the exact solution, the BCS-pair order parameters, with respect to two independent sub-systems respectively, are introduced to characterize the phase diagram by its value and nonanalytic behavior at phase boundaries. We find that the ground state are two-fluid condensate, with two different BCS pairing modes. In addition, such a decomposition holds true not only for ground state but also the whole eigen space. Then it should result in the decomposition of dynamics. Based on the exact results on the dynamics obtained in the previous workSYB_PRB , we study the dynamics of the SSH Kitaev model. The primary aim of this study was to obtain pertinent information regarding the nonequilibrium state. In this study, the starting point is the vacuum state, a trivial initial state. We investigate the evolved states under the postquench Hamiltonian with parameters covering the entire region and find two types of order parameters that are different but can help determine the quantum phase diagram. Finally, the robustness of such a decomposition to the perturbation of nonzero chemical potential is also investigated. Analytical analysis and numerical simulations show that the real-space decomposition for the ground state still holds true approximately in many regions in the presence of finite chemical potential.

This paper is organized as follows. In Sec. II, we present a generalized model and introduce a real-space Bogoliubov transformation to decompose the model Hamiltonian. In Secs. III, we apply the transformation on a SSH Kitaev model and deduce the phase diagram from decomposed sub-Hamiltonians. In Secs. IV and V, two types of BCS order parameters are introduced. The explicit expressions for both ground state and non-equilibrium state are obtained, respectively. Sec. VI devotes to the investigation for the case with nonzero potential. Finally, we provide a summary and discussion in Section VII.

We consider the following fermionic Hamiltonian on a lattice of length $2N$

$$H=\sum\limits_{j=1}^{2N}[J_{j}c_{j}^{{\dagger}}c_{j+1}+\Delta_{j}c_{j}^{{\dagger}}c_{j+1}^{{\dagger}}+\mathrm{H.c.}+\mu_{j}\left(2n_{j}-1\right)],$$

(1)

where $c_{j}^{{\dagger}}$ $(c_{j})$ is a fermionic creation (annihilation) operator on site $j$, $n_{j}=c_{j}^{{\dagger}}c_{j}$, $J_{j}$ the tunneling rate across the dimer $\left(j,j+1\right)$, $\mu_{j}$ the on-site chemical potential, and $\Delta_{j}$ the strength of the $p$-wave pair creation (annihilation) on the dimer $\left(j,j+1\right)$. $c_{2N+1}=c_{1}$ is defined for periodic boundary condition. The Hamiltonian in Eq. (1) has a rich phase diagram that describes a spin-polarized $p$-wave superconductor in one dimension when the parameters $\left(J_{j},\Delta_{j},\mu_{j}\right)$ are taken in a uniform fashion Ryohei Wakatsuki . This system has a topological phase in which a zero- energy Majorana mode is located at each end of a long chain by taking the open boundary condition $c_{2N+1}=0$. It is the fermionized version of the well-known one-dimensional transverse-field Ising model Pfeuty , which is one of the simplest solvable models when the translational symmetry is imposed that exhibits quantum criticality and phase transition with spontaneous symmetry breaking SachdevBook ; ZG . In addition, several studies have been conducted with a focus on long-range Kitaev chains DV1 ; DV2 ; OV ; LL ; UB . Although system we concerned in this work contains only the nearest neighbor hopping and pairing terms, this method be applied to the system involving long-range terms under certain conditions.

We introduce a linear unitary transformation

$$\left(\begin{array}[]{c}\alpha_{j}^{{\dagger}}\\ \alpha_{j}\\ \beta_{j}^{{\dagger}}\\ \beta_{j}\end{array}\right)=T\left(\begin{array}[]{c}c_{2j-1}^{{\dagger}}\\ c_{2j-1}\\ c_{2j}^{{\dagger}}\\ c_{2j}\end{array}\right)=T\left(\begin{array}[]{c}A_{j}^{{\dagger}}\\ A_{j}\\ B_{j}^{{\dagger}}\\ B_{j}\end{array}\right),$$

(2)

with the matrix

$$T=\left(T^{{\dagger}}\right)^{-1}=\frac{1}{2}\left(\begin{array}[]{cccc}1&1&1&-1\\ 1&1&-1&1\\ -i&i&i&i\\ -i&i&-i&-i\end{array}\right),$$

(3)

which ensures that the $\alpha_{j}$ and $\beta_{j}$ are still fermionic operators, i.e.,

	$\displaystyle\left\{\alpha_{j},\alpha_{j^{\prime}}^{{\dagger}}\right\}$	$\displaystyle=$	$\displaystyle\left\{\alpha_{j},\beta_{j^{\prime}}^{{\dagger}}\right\}=\delta_{j,j^{\prime}},\left\{\alpha_{j},\beta_{j^{\prime}}^{{\dagger}}\right\}=0,$
	$\displaystyle\left\{\alpha_{j},\alpha_{j^{\prime}}\right\}$	$\displaystyle=$	$\displaystyle\left\{\beta_{j},\beta_{j^{\prime}}\right\}=\left\{\alpha_{j},\beta_{j^{\prime}}\right\}=0,$		(4)

and any transformed Hamiltonian maintains its original physics. It is similar to the Bogoliubov transformation if regarding the parity of the site number ($2j$ or $2j-1$) as spin degree of freedom. Applying the transformation (2), the Hamiltonian $H$ can be expressed in the form

$$H=H_{A}+H_{B}+H_{AB},$$

(5)

where

	$\displaystyle H_{A}=\frac{1}{2}\sum\limits_{j=1}^{N}[\left(J_{2j}+\Delta_{2j}\right)\left(\alpha_{j}^{{\dagger}}\alpha_{j+1}+\alpha_{j}^{{\dagger}}\alpha_{j+1}^{{\dagger}}\right)$
	$\displaystyle+\mathrm{H.c.}+\left(\Delta_{2j-1}-J_{2j-1}\right)\left(2\alpha_{j}^{{\dagger}}\alpha_{j}-1\right)],$		(6)

and

	$\displaystyle H_{B}=\frac{1}{2}\sum\limits_{j=1}^{N}[\left(\Delta_{2j}-J_{2j}\right)\left(\beta_{j}^{{\dagger}}\beta_{j+1}+\beta_{j}^{{\dagger}}\beta_{j+1}^{{\dagger}}\right)$
	$\displaystyle+\mathrm{H.c.}+\left(J_{2j-1}+\Delta_{2j-1}\right)\left(2\beta_{j}^{{\dagger}}\beta_{j}-1\right)],$		(7)

denote two Kitaev chain of length $N$ and commute with each other

$$\left[H_{A},H_{B}\right]=0.$$

(8)

The third term is $\left\{\mu_{l}\right\}$ dependent, containing the intra- and inter-chain interactions

	$\displaystyle H_{AB}$	$\displaystyle=$	$\displaystyle i\sum_{j=1}^{N}[\left(\mu_{2j-1}+\mu_{2j}\right)\beta_{j}^{{\dagger}}\alpha_{j}^{{\dagger}}$		(9)
			$\displaystyle+\left(\mu_{2j-1}-\mu_{2j}\right)\beta_{j}^{{\dagger}}\alpha_{j}]+\mathrm{H.c.}$

Obviously, when consider the case with zero chemical potential, $\mu_{l}=0$, the Hamiltonian is exactly decomposed into two independent ones. Accordingly we have

$$H\left|\psi\right\rangle=H\left|\psi_{A}\right\rangle\left|\psi_{B}\right\rangle=\left(\varepsilon_{A}+\varepsilon_{B}\right)\left|\psi\right\rangle,$$

(10)

with

$$H_{A}\left|\psi_{A}\right\rangle=\varepsilon_{A}\left|\psi_{A}\right\rangle,H_{B}\left|\psi_{B}\right\rangle=\varepsilon_{B}\left|\psi_{B}\right\rangle.$$

(11)

For arbitrary operator functions $F_{A}(\alpha_{1},...,$ $\alpha_{l},...,\alpha_{N};$ $\alpha_{1}^{{\dagger}},...,$ $\alpha_{j}^{{\dagger}},...,\alpha_{N}^{{\dagger}})$ and $F_{B}(\beta_{1},...,$ $\beta_{l},...,$ $\beta_{N};$ $\beta_{1}^{{\dagger}},...,$ $\beta_{j}^{{\dagger}},...,\beta_{N}^{{\dagger}})$, we always have

	$\displaystyle\left\langle\psi\right|F_{A}\left|\psi\right\rangle$	$\displaystyle=$	$\displaystyle\left\langle\psi_{A}\right|F_{A}\left|\psi_{A}\right\rangle,$
	$\displaystyle\left\langle\psi\right|F_{B}\left|\psi\right\rangle$	$\displaystyle=$	$\displaystyle\left\langle\psi_{B}\right|F_{B}\left|\psi_{B}\right\rangle.$		(12)

Then any features relate to $F_{A}$($F_{B}$) in chain $A$($B$) should emerge in the original system. The the following we will demonstrate this in a concrete example.

Now we consider a SSH Kitaev model with the Hamiltonian

	$\displaystyle H=\sum_{j=1}^{N}[\lambda c_{2j}^{{\dagger}}c_{2j+1}+\left(1-\lambda\right)c_{2j-1}^{{\dagger}}c_{2j}+\lambda\Delta c_{2j-1}c_{2j}$
	$\displaystyle+\left(1-\lambda\right)\Delta c_{2j}c_{2j+1}+2\mu\left(n_{2j-1}+n_{2j}-1\right)]+\mathrm{H.c.},$		(13)

in which, both the hopping and pairing strengths are assigned alternatively. It has been well studied in the refs. Ryohei Wakatsuki . In this work, we focus on the connection between this $2N$-site system and two $N$-site sub-systems. When taking $\mu=0$, the original Hamiltonian can be written in the form $H=H_{A}+H_{B}$ with

	$\displaystyle H_{\sigma}$	$\displaystyle=$	$\displaystyle\sum_{j=1}^{N}[J_{\sigma}\left(d_{j,\sigma}^{{\dagger}}d_{j+1,\sigma}^{{\dagger}}+d_{j,\sigma}^{{\dagger}}d_{j+1,\sigma}\right)+\text{H.c.}$		(14)
			$\displaystyle+\mu_{\sigma}\left(1-2d_{j,\sigma}^{{\dagger}}d_{j,\sigma}\right),$

where $\sigma=A$ or $B$, $d_{j,A}=\alpha_{j}$, $d_{j,B}=\beta_{j}$, and the corresponding parameters are

	$\displaystyle J_{A}$	$\displaystyle=$	$\displaystyle\left[\lambda-\left(1-\lambda\right)\Delta\right]/2,$		(15)
	$\displaystyle\mu_{A}$	$\displaystyle=$	$\displaystyle\left[\lambda\Delta+\left(1-\lambda\right)\right]/2,$		(16)
	$\displaystyle J_{B}$	$\displaystyle=$	$\displaystyle-\left[\lambda+\left(1-\lambda\right)\Delta\right]/2,$		(17)
	$\displaystyle\mu_{B}$	$\displaystyle=$	$\displaystyle\left[\lambda\Delta-\left(1-\lambda\right)\right]/2.$		(18)

We note that both $H_{A}$ and $H_{B}$ describe the same system but with different parameters. Although the chemical potential is zero in $H$, it is nonzero for sub-Hamiltonians. The phase diagram for each sub-Hamiltonian $H_{\sigma}$ is well known, and then can be used to obtain the phase diagram of $H$. In fact, the phase boundary of the ground state of $H_{\sigma}$ are two points $\mu_{\sigma}/J_{\sigma}=\pm 1$. This maps to two curves

$$\Delta-2\lambda+1=0,$$

(19)

and

$$2\lambda\Delta-\Delta+1=0,$$

(20)

in the $\lambda\Delta$-plane for $H_{A}$, while two curves

$$\Delta+2\lambda-1=0,$$

(21)

and

$$2\lambda\Delta-\Delta-1=0,$$

(22)

for $H_{B}$. According to our above analysis, all the four curves constitute the phase diagram of $H$. In Fig. 1, we schematically illustrate the phase diagrams of $H_{A}$, $H_{B}$ and $H$, respectively. We can see that all the information in the ground state of $H$, including the gapless lines and winding numbers can be obtained directly by that of $H_{A}$ and $H_{B}$, as a demonstration of benefits from the real-space decomposition.

It has been shown that the pairing order parameter can be utilized to characterize the phase diagram SYB_PRB and the long-range order MES_PRB of the ground state of a simple Kitaev chain. In this section, we focus on the similar investigation in this aspect for the present model based on the real-space decomposition method. We first briefly review the obtained conclusion for the simple Kitaev chain. We take $H_{A}$ as an example. We introduce the BCS-pairing operator

$$\widehat{O}_{A,k}=i\left(\alpha_{-k}\alpha_{k}-\alpha_{k}^{{\dagger}}\alpha_{-k}^{{\dagger}}\right),$$

(23)

to characterize pairing channels in ${k}$ space. Here the Fourier transformation

	$\displaystyle\alpha_{k}$	$\displaystyle=$	$\displaystyle\frac{1}{\sqrt{N}}\sum\limits_{j}e^{-ikj}\alpha_{j}$		(24)
		$\displaystyle=$	$\displaystyle\frac{1}{2}\left(A_{k}+B_{k}+A_{-k}^{{\dagger}}-B_{-k}^{{\dagger}}\right),$

is applied with

$$\left(\begin{array}[]{c}A_{k}\\ B_{k}\end{array}\right)=\frac{1}{\sqrt{N}}\sum_{k}e^{-ikj}\left(\begin{array}[]{c}A_{j}\\ B_{j}\end{array}\right).$$

(25)

For a given state $\left|\psi\right\rangle$, the quantity $|\left\langle\psi\right|\widehat{O}_{A,k}\left|\psi\right\rangle|$ measures the rate of transition for a BCS pair at the $k$ channel. For the ground state $\left|\text{{G}}_{A}\right\rangle$ of $H_{A}$ the pairing order parameter of the ground state is expressed as

$$O_{A,\mathrm{g}}=\frac{1}{N}\sum_{\pi>k>0}\left|\left\langle\text{{G}}_{A}\right|\widehat{O}_{A,k}\left|\text{{G}}_{A}\right\rangle\right|.$$

(26)

In the large $N$ limit, it can be expressed explicitly as

	$\displaystyle O_{A,\mathrm{g}}$	$\displaystyle=$	$\displaystyle\frac{1}{\pi}\int_{0}^{\pi}\frac{\sin k}{2\sqrt{\left(\mu_{A}/J_{A}-\cos k\right)^{2}+\sin^{2}k}}\mathrm{d}k$		(29)
		$\displaystyle=$	$\displaystyle\frac{1}{\pi}\left\{\begin{array}[]{cc}1,&\left|\mu_{A}/J_{A}\right|\leqslant 1\\ J_{A}/\mu_{A},&\left|\mu_{A}/J_{A}\right|>1\end{array}\right..$

Obviously, there exist non-analytic points at $\left|\mu_{A}/J_{A}\right|=1$, as the signatures of quantum phase boundary. The same results hold for $H_{B}$.

These conclusions can be directly applied on the ground state $\left|\text{{G}}\right\rangle$ of system $H$, mapping the non-analytic points $\left|\mu_{A}/J_{A}\right|=\left|\mu_{B}/J_{B}\right|=1$ on the $\lambda\Delta$-plane. In fact, introducing two order parameters

$$O_{A,\mathrm{g}}=\frac{1}{N}\sum_{\pi>k>0}\left|\left\langle\text{{G}}\right|i\left(\alpha_{-k}\alpha_{k}-\alpha_{k}^{{\dagger}}\alpha_{-k}^{{\dagger}}\right)\left|\text{{G}}\right\rangle\right|,$$

(30)

and

$$O_{B,\mathrm{g}}=\frac{1}{N}\sum_{\pi>k>0}\left|\left\langle\text{{G}}\right|i\left(\beta_{-k}\beta_{k}-\beta_{k}^{{\dagger}}\beta_{-k}^{{\dagger}}\right)\left|\text{{G}}\right\rangle\right|,$$

(31)

we directly have the explicit expressions of $O_{A,\mathrm{g}}$ and $O_{B,\mathrm{g}}$ in the $\lambda\Delta$-plane

$$O_{A,\mathrm{g}}=\frac{1}{\pi}\left\{\begin{array}[]{cc}1,&\text{regions I, II, and III}\\ \frac{\lambda-\left(1-\lambda\right)\Delta}{\lambda\Delta+\left(1-\lambda\right)},&\text{otherwise}\end{array}\right.,$$

(32)

where region I is $\Delta\in(-\frac{1}{2\lambda-1},\infty)$ for $0<\lambda<0.5$, region II is $\Delta\in(-\frac{1}{2\lambda-1},2\lambda-1)$ for $0.5<\lambda<1$, and region III is $\Delta\in(-\infty,2\lambda-1)$ for $0<\lambda<0.5$, while

$$O_{B,\mathrm{g}}=\frac{1}{\pi}\left\{\begin{array}[]{cc}1,&\text{regions I, II, and III}\\ -\frac{\lambda+\left(1-\lambda\right)\Delta}{\lambda\Delta-\left(1-\lambda\right)},&\text{otherwise}\end{array}\right.,$$

(33)

where region I is $\Delta\in(-\infty,\frac{1}{2\lambda-1})$ for $0<\lambda<0.5$, region II is $\Delta\in(-2\lambda+1,\frac{1}{2\lambda-1})$ for $0.5<\lambda<1$, and region III is $\Delta\in(-2\lambda+1,\infty)$ for $0<\lambda<0.5$. In the representation of the original Hamiltonian, the physics of the two order parameters are clear, based on the expressions

	$\displaystyle O_{A,\mathrm{g}}+O_{B,\mathrm{g}}$	$\displaystyle=$	$\displaystyle\frac{1}{N}\sum_{\pi>k>0}\left|\left\langle\text{{G}}\right|\Lambda_{k}\left|\text{{G}}\right\rangle\right|,$		(34)
	$\displaystyle O_{A,\mathrm{g}}-O_{B,\mathrm{g}}$	$\displaystyle=$	$\displaystyle\frac{1}{N}\sum_{\pi>k>0}\left|\left\langle\text{{G}}\right|\Sigma_{k}\left|\text{{G}}\right\rangle\right|,$		(35)

with the pairing and current operators

	$\displaystyle\Lambda_{k}$	$\displaystyle=$	$\displaystyle i\left(A_{-k}B_{k}+B_{-k}A_{k}-\text{{H.c.}}\right),$		(36)
	$\displaystyle\Sigma_{k}$	$\displaystyle=$	$\displaystyle i\left(A_{k}^{{\dagger}}B_{k}+B_{-k}A_{-k}^{{\dagger}}-\text{{H.c.}}\right).$		(37)

Obviously, $O_{A,\mathrm{g}}+O_{B,\mathrm{g}}$ measures the BCS pair for fermions from two sub-lattices $A$ and $B$, while $O_{A,\mathrm{g}}-O_{B,\mathrm{g}}$ the current across two sub-lattices. In Fig. 2, we plot the functions, $O_{A,\mathrm{g}}$, $O_{B,\mathrm{g}}$ and $O_{A,\mathrm{g}}\pm O_{B,\mathrm{g}}$, respectively. We can see that quantities $O_{A,\mathrm{g}}$ and $O_{B,\mathrm{g}}$ reflect partial phase boundary, while both two types of order parameters $O_{A,\mathrm{g}}\pm O_{B,\mathrm{g}}$ can identify the entire phase diagram.

Recently, advancements in atomic physics, quantum optics, and nanoscience have allowed the development of artificial systems with high accuracy Jochim ; Greiner . The study of nonequilibrium many-body dynamics presents an alternative approach for accessing a new exotic quantum state with an energy level considerably different from that of the ground state Choi ; Else ; Khemani ; Lindner ; Kaneko ; Tindall ; YXMPRA ; ZXZPRB2 ; TK ; JT1 ; JT2 . Several works show that the quenching dynamics governed by the post-quench Hamiltonian is intimately related to its ground state SYB_PRB ; MH ; LZ . In this section, we turn to the topic of nonequilibrium phenomena in the present Hamiltonian. For a uniform Kitaev chain, it has been shown that the pairing order parameter for a nonequilibrium state obtained through time evolution from an initially prepared vacuum state still help determine the phase diagram SYB_PRB . For the present model, we will investigate this issue based on the obtained results for the sub-systems $H_{A}$ and $H_{B}$.

Similarly, we still give a brief review the obtained conclusion for a simple Kitaev chain, such as chain $H_{A}$. It has been shown that the order parameter for a nonequilibrium state is defined as SYB_PRB

$$\overline{O}_{A}=\lim_{T\rightarrow\infty}\frac{1}{T}\int_{0}^{T}\frac{1}{N}\sum_{k}\left|\left\langle\psi(t)\right|\widehat{O}_{A,k}\left|\psi(t)\right\rangle\right|\mathrm{d}t,$$

(38)

where the time evolution $\left|\psi_{A}(t)\right\rangle=e^{-iH_{A}t}\left|\psi_{A}(0)\right\rangle$ for a particular initial state $\left|\psi(0)\right\rangle$, satisfying

$$\alpha_{j}\left|\psi_{A}(0)\right\rangle=0.$$

(39)

State $\left|\psi(0)\right\rangle$ is essentially an empty state for a set of fermions $\left\{\alpha_{j}\right\}$ in real space. In the large $N$ limit, it can be expressed explicitly as

	$\displaystyle\overline{O}_{A}$	$\displaystyle=$	$\displaystyle\frac{1}{2\pi}\int_{0}^{\pi}\left|\frac{\left(\cos k-\mu_{A}/J_{A}\right)\sin k}{\left(\cos k-\mu_{A}/J_{A}\right)^{2}+\sin^{2}k}\right|\mathrm{d}k$		(42)
		$\displaystyle=$	$\displaystyle\frac{1}{2\pi}\left\{\begin{array}[]{cc}1,&\left|\mu_{A}/J_{A}\right|\leqslant 1\\ \Lambda,&\left|\mu_{A}/J_{A}\right|>1\end{array}\right.,$

where

$$\Lambda=\left|\frac{J_{A}}{\mu_{A}}+\frac{1}{2}\left(1-\frac{J_{A}^{2}}{\mu_{A}^{2}}\right)\ln\left|\frac{\mu_{A}+J_{A}}{\mu_{A}-J_{A}}\right|\right|.$$

(43)

The non-analytic points at $\left|\mu_{A}/J_{A}\right|=1$ are the evidently signatures of quantum phase boundary. The same results hold for $H_{B}$.

These conclusions can be directly applied to the dynamics of system $H$, mapping the non-analytic points $\left|\mu_{A}/J_{A}\right|=\left|\mu_{B}/J_{B}\right|=1$ on the $\lambda\Delta$-plane. Accordingly, the nonequilibrium order parameters are defined for the time evolution $\left|\psi(t)\right\rangle=e^{-iHt}\left|\psi(0)\right\rangle$ $=e^{-iH_{A}t}e^{-iH_{B}t}\left|\psi(0)\right\rangle$ of a particular initial state $\left|\psi(0)\right\rangle$. To utilize the conclusion for the two sub-systems, the initial state should be chosen as the vacuum state $\left|\psi(0)\right\rangle=\left|\text{{Vac}}\right\rangle$ for both sets of fermion operators $\left\{\alpha_{j}\right\}$ and $\left\{\beta_{j}\right\}$, i.e.,

$$\alpha_{j}\left|\text{{Vac}}\right\rangle=\beta_{j}\left|\text{{Vac}}\right\rangle=0,$$

rather than the vacuum state of operator $\left\{c_{j}\right\}$. The vacuum state can be constructed as the form

$$\left|\text{{Vac}}\right\rangle=\prod_{j=1}^{N}\beta_{j}\alpha_{j}\left|0\right\rangle=\prod_{j=1}^{N}\frac{i}{\sqrt{2}}\left(1+c_{2j-1}^{{\dagger}}c_{2j}^{{\dagger}}\right)\left|0\right\rangle,$$

(44)

which is essentially the ground state of the Hamiltonian $H$ at the special point with $\lambda=1$ and $\Delta\ll-1$. Then such a special vacuum state can be prepared in this way in the experiment.

By replacing the parameters $\left\{\mu_{A}/J_{A},\mu_{B}/J_{B}\right\}$ by $\left\{\Delta,\lambda\right\}$, we obtain the explicit expressions of $\overline{O}_{A}$ and $\overline{O}_{A}$ in the $\lambda\Delta$-plane

$$\overline{O}_{A}=\frac{1}{2\pi}\left\{\begin{array}[]{cc}1,&\text{regions I, II, and III}\\ \Omega_{A},&\text{otherwise}\end{array}\right.$$

(45)

where region I is $\Delta\in(-\frac{1}{2\lambda-1}$, $\infty)$ for $0<\lambda<0.5$, region II is $\Delta\in(-\frac{1}{2\lambda-1}$, $2\lambda-1)$ for $0.5<\lambda<1$, and region III is $\Delta\in(-\infty$, $2\lambda-1)$ for $0<\lambda<0.5$,

$$\overline{O}_{B}=\frac{1}{2\pi}\left\{\begin{array}[]{cc}1,&\text{regions I, II, and III}\\ \Omega_{B},&\text{otherwise}\end{array}\right.$$

(46)

where region I is $\Delta\in(-\infty$, $\frac{1}{2\lambda-1})$ for $0<\lambda<0.5$, region II is $\Delta\in(-2\lambda+1$, $\frac{1}{2\lambda-1})$ for $0.5<\lambda<1$, and region III is $\Delta\in(-2\lambda+1$, $\infty)$ for $0<\lambda<0.5$. Here two factors can be expressed explicitly as

$$\Omega_{A}=\left|\frac{\lambda-\left(1-\lambda\right)\Delta}{\lambda\Delta+\left(1-\lambda\right)}+\frac{1}{2}\left[1-\left(\frac{\lambda-\left(1-\lambda\right)\Delta}{\lambda\Delta+\left(1-\lambda\right)}\right)^{2}\right]\ln\left|\frac{2\lambda\Delta+1-\Delta}{1-2\lambda+\Delta}\right|\right|,$$

(47)

and

$$\Omega_{B}=\left|-\frac{\lambda+\left(1-\lambda\right)\Delta}{\lambda\Delta-\left(1-\lambda\right)}+\frac{1}{2}\left[1-\left(\frac{\lambda+\left(1-\lambda\right)\Delta}{\lambda\Delta-\left(1-\lambda\right)}\right)^{2}\right]\ln\left|\frac{2\lambda\Delta-1-\Delta}{-1+2\lambda+\Delta}\right|\right|.$$

(48)

In Fig. 3, we plot the functions, $\overline{O}_{A}$, $\overline{O}_{B}$ and $\overline{O}_{A}\pm\overline{O}_{B}$, respectively. We can see that quantities $\overline{O}_{A}$ and $\overline{O}_{B}$ reflect partial phase boundary, while both two types of order parameters $\overline{O}_{A}\pm\overline{O}_{B}$ can identify the entire phase diagram. Methodologically, it also demonstrates the benefits of the real-space decomposition.

The real-space decomposition is exact when the chemical potential is zero. In this section, we investigate the influence of nonzero $\mu$ on our above conclusions. The decomposition we proposed is no longer valid, in the presence of nonzero $\mu$. However, we will show that our conclusion for the ground state holds approximately in most of regions for small values of $\mu$. In principle, when small $\mu$ switches on, one can treat $H_{AB}$ as a perturbation term. The perturbation theory tells us the effect of the perturbation term depends strongly on the eigenstates of $H_{A}+H_{B}$. In the limit case, the energy gaps for the groundstate of both $H_{A}+H_{B}$ are sufficiently large comparing to the value of $\mu$, $H_{AB}$ should have no effect on the groundstates $\left|\text{{G}}_{A}\right\rangle\left|\text{{G}}_{B}\right\rangle$ of $H_{A}+H_{B}$. On the other hand, at the phase boundary, the gapless point, the ground states $\left|\text{{G}}_{A}\right\rangle\left|\text{{G}}_{B}\right\rangle$ become (or quasi-) degenerate. The term $H_{AB}$ may hybridize them and the low-lying excited states, resulting in entangled state, which is deviated from the product state $\left|\text{{G}}_{A}\right\rangle\left|\text{{G}}_{B}\right\rangle$.

Based on the Fourier transformations for two sub-lattices in Eq. (25), the Hamiltonian with periodic boundary condition can be block diagonalized by this transformation due to its translational symmetry, i.e.,

$$H=\sum_{k\in\left[0,\pi\right]}H_{k}=H_{0}+H_{\pi}+\sum_{k\in(0,\pi)}\psi_{k}^{\dagger}h_{k}\psi_{k},$$

(49)

satisfying $\left[H_{k},H_{k^{\prime}}\right]=0$, where the operator vector $\psi_{k}^{\dagger}=\left(A_{k}^{\dagger},B_{k}^{\dagger},A_{-k},B_{-k}\right)$, and the core matrix is expressed explicitly as

$$h_{k}=\left(\begin{array}[]{cccc}\mu&z&0&w\\ z^{\ast}&\mu&-w^{\ast}&0\\ 0&-w&-\mu&-z\\ w^{\ast}&0&-z^{\ast}&-\mu\end{array}\right),$$

(50)

where $z=\lambda e^{ik}$ $+\left(1-\lambda\right)$ and $w=\Delta\left[\lambda\text{ }-\left(1-\lambda\right)e^{ik}\right]$. Here, $H_{0}$ and $H_{\pi}$ have the form

	$\displaystyle H_{0}$	$\displaystyle=$	$\displaystyle 2JA_{0}^{\dagger}B_{0}+2JB_{0}^{\dagger}A_{0}+2\mu\left(A_{0}^{{\dagger}}A_{0}-B_{0}B_{0}^{{\dagger}}\right)$		(51)
			$\displaystyle+\left(\Delta_{a}-\Delta_{b}\right)\left(A_{0}^{\dagger}B_{0}^{\dagger}+A_{0}B_{0}\right),$
	$\displaystyle H_{\pi}$	$\displaystyle=$	$\displaystyle 2\mu\left(A_{\pi}^{{\dagger}}A_{\pi}-B_{\pi}B_{\pi}^{{\dagger}}\right)$		(52)
			$\displaystyle-\left(\Delta_{a}+\Delta_{b}\right)\left(A_{\pi}^{\dagger}B_{\pi}^{\dagger}+B_{\pi}A_{\pi}\right),$

and $H_{\pi}$ vanishes when odd $N$ is taken. To demonstrate the above analysis, one can rewrite the matrix $h_{k}$ in the form

$$h_{k}=h_{k}^{0}+\mu\Gamma^{z}$$

(53)

with

$$h_{k}^{0}=\left(\begin{array}[]{cccc}0&z&0&w\\ z^{\ast}&0&-w^{\ast}&0\\ 0&-w&0&-z\\ w^{\ast}&0&-z^{\ast}&0\end{array}\right),$$

(54)

and a $k$-independent matrix

$$\Gamma^{z}=\left(\begin{array}[]{cccc}1&0&0&0\\ 0&1&0&0\\ 0&0&-1&0\\ 0&0&0&-1\end{array}\right).$$

(55)

The eigen values $\varepsilon_{\sigma\rho}^{k}$ and vectors $\phi_{\sigma\rho}^{k}$ of $h_{k}^{0}$ can be obtained as

$$\left(\begin{array}[]{c}\left(\phi_{++}^{k}\right)^{T}\\ \left(\phi_{+-}^{k}\right)^{T}\\ \left(\phi_{-+}^{k}\right)^{T}\\ \left(\phi_{--}^{k}\right)^{T}\end{array}\right)=\frac{1}{2}\left(\begin{array}[]{cccc}\eta_{k}&1&-\eta_{k}&1\\ \eta_{k}&-1&\eta_{k}&1\\ \eta_{k}&-1&-\eta_{k}&-1\\ \eta_{k}&1&\eta_{k}&-1\end{array}\right),$$

(56)

and

$$\varepsilon_{\sigma\rho}^{k}=\sigma\left|w+\rho z\right|$$

(57)

satisfying $h_{k}^{0}\phi_{\sigma\rho}^{k}=\varepsilon_{\sigma\rho}^{k}\phi_{\sigma\rho}^{k}$, where $\eta_{k}=\frac{w+z}{\left|w+z\right|}$, with the indices $\sigma,\rho=\pm$. It is easy to check that the $h_{k}^{\mu}$ acts as a flip operator

$$\Gamma^{z}\phi_{\sigma\rho}^{k}=\phi_{\overline{\sigma}\overline{\rho}}^{k},$$

(58)

with the labels $\left(\overline{\sigma},\overline{\rho}\right)=\left(-\sigma,-\rho\right)$. It directly results in the transition matrix element between energy bands

$$\left(\phi_{\sigma\rho}^{k}\right)^{{\dagger}}\Gamma^{z}\phi_{\sigma^{\prime}\rho^{\prime}}^{k}=\delta_{\sigma^{\prime}\overline{\sigma}}\delta_{\rho^{\prime}\overline{\rho}}.$$

(59)

This indicates that when the energy gap is sufficiently large compared to the values of the chemical potential,

$$\left|\mu\right|\ll\left|\varepsilon_{\sigma\rho}^{k}-\varepsilon_{\overline{\sigma}\overline{\rho}}^{k}\right|,$$

(60)

the term $\mu\Gamma^{z}$ should not affect the eigen vectors of $h_{k}^{0}$.

We introduce the quantity

$$f\left(\lambda,\Delta\right)=\left|\langle\text{{G}}(\mu)\left|\text{{G}}_{A}\right\rangle\left|\text{{G}}_{B}\right\rangle\right|^{2}=\left|\langle\text{{G}}(\mu)\left|\text{{G}}(0)\right\rangle\right|^{2},$$

(61)

to quantitatively measure the effect of $\mu$ on the groundstate $\left|\text{{G}}\right\rangle$ of $H$. Our above analysis predicts that the gapless lines of both $H_{A}$ and $H_{B}$ result in the valley of $f\left(\lambda,\Delta\right)$. Numerical simulations for $f\left(\lambda,\Delta\right)$ in finite size system can be performed by exact diagonalization. It essentially relates to two ground states of $H$ with zero and nonzero $\mu$, respectively, while $\left|\text{{G}}(0)\right\rangle=\left|\text{{G}}_{A}\right\rangle\left|\text{{G}}_{B}\right\rangle$ is the ground state for zero $\mu$. For given parameters, $f\left(\lambda,\Delta\right)$ can be obtained by diagonalization of matrix $h_{k}$ and the eigenstates of $H_{0}+H_{\pi}$.

Fig. 4 shows the fidelity of $f\left(\lambda,\Delta\right)$ for a given finite system as a function of $\left(\lambda,\Delta\right)$ with various values of $\mu$. In the case of zero $\mu$, the fidelity is unitary everywhere as expected. In plot for small $\mu$, there are evidently sudden drops of the value of fidelity at the phase boundaries and the drops become sharper and sharper as $\mu$ increases. These behaviors can be ascribed to a dramatic change of the ground state of the system around the phase boundaries. One can also find out that the ground state remains unchanged within many regions. It indicates the real-space decomposition for the ground state still holds true approximately in many regions in the presence of finite $\mu$.

In summary, we have extended the Bogoliubov transformation for spinless fermions in $k$-space to the one in real space, which has been shown to be a tool for the spatial decomposition of a class of Kitaev chains. It provides a way to get insight into a complicated system from that of two decoupled sub-systems. A systematic investigation of a SSH Kitaev model is performed, including the ground state property and nonequilibrium behavior of quenching dynamics. In addition, we also studied the approximate decomposition in the case of nonzero chemical potential. Analytical analysis and numerical simulation show that the real-space decomposability is maintained in most regions of the phase diagram when the chemical potential is small. Our findings not only contribute to the methodology for solving the many-body problem, but also reveal the underlying mechanism of the features of the SSH Kitaev chain.

We acknowledge the support of NSFC (Grants No. 12374461).