Continuous phase transition induced by non-Hermiticity in the quantum contact process model

Non-Hermitian physics currently has attracted considerable attentions [1, 2], which greatly extend the knowledge about quantum system. For example, the complex energy coalesce at exceptional point(EP)[2, 3, 4, 5]; The wave function localize at the edge of quantum system, which is now named as non-Hermitian skin effect[3, 4, 6, 7]. Actually, quantum systems almost are subjected to environment. It is more appropriate to consider the quantum system to be non-Hermitian(NH). Recently, it was found that driven-dissipative non-Hermitian system can provide one method to control quantum system, for instance, directional amplification [8, 9, 10, 11]. Non-Hermiticity can impose great influence on the properties of the topological non-trivial system. It needs generalized Brillouin zone to establish a generalized Bloch band theory [4, 12, 13]. Non-Hermiticity can also induce phase transition without gap[14] and non-Hermitian many-body localization [15]. Non-Hermitian quantum system needs more effort to be studied, understood and utilized.

Though we have understood several interesting property of non-Hermitian quantum system which based on single-particle picture, it is an important issue to understand that how non-Hermiticity affect many-body system, especially the phase transition of many-body system. The theory of phase transition provide people profound understanding to the status of the matter. When non-Hermiticity comes into quantum many-body system, whether non-Hermiticity can induce phase transition, and if the phase transition exists how to classify this type of phase transition. These questions remain being elusive. Moreover, what the relation is between the phase transition induced by non-Hermiticity and their Hermitian counterpart. As the primary motivation to answer above questions, we are going to explore the phase transition induced by non-Hermiticity in quantum contact process model.

In this work, we consider the open system whose dynamics governed by Lindbladian master equation. We use effective Hamiltonian to describe the property of non-Hermitian quantum system. Then we choose quantum contact process(QCP) model as research object, which breaks $U(1)$ symmetry. We combine analytical method for special case of spin number $L=2$ and exact diagonalization(ED) for general $L$ to compute energy spectrum, the order parameter $M^{x}$, susceptibility $\chi$, correlation, and entanglement entropy. The results of this work show that there is continuous phase transition within the non-Hermitian QCP model. We extract the critical exponents, $\beta$ of order parameter, $\gamma$ of susceptibility. We reveal the transition from quasi long-rang order to short-range order and entanglement properties near critical point. Finally, we compare the results of non-Hermitian QCP model with Hermitian counterpart, which shows the phase transition in non-Hermitian QCP model have no correspondence in Hermitian case. By looking at the critical exponents, the results in our work may indicate the phase transition belongs a new universality class. Moreover, even for finite size system, we notice non-Hermiticity endow many-body system with different singular behaviour from classical phase transition. Our results can extend the knowledge about the phase transition of non-Hermitian quantum matter.

For the open quantum many-body systems, its dynamics usually can be determined by Lindbladian master equation [16] as

$$\dot{\rho}(t)=-i[\hat{H}_{0},\rho(t)]+\sum_{k}^{L}{{\cal D}_{k}[\rho(t)]}.$$

(1)

where $\hat{H}_{0}$ is Hamiltonian governing the coherent evolution of the system. The dissipation superoperator ${\cal D}_{\alpha}[\rho(t)]=\hat{L}^{k}\rho(t)\hat{L}^{k,\dagger}-\frac{1}{2}\{\hat{L}^{k^{\dagger}}\hat{L}^{k},\rho(t)\}$ is defined by the following Lindbladian jump operator, $\hat{L}^{k}=\sqrt{\Gamma}\sigma_{-}^{k}$ which impose local spin decaying channels. In the absence of quantum jumps during a time interval, the term $L^{k}\rho(t)L^{k\dagger}$, which describes the state transition, can be omitted[1, 17]. This corresponds to postselection [17, 18]in the experiment, which target on succeed measurement results. In this case, the quantum system is considered to undergo a non-unitary time evolution governed by an effective Hamiltonian which is derived from the master equation(1) as follows

$$\hat{H}_{\rm eff}=\hat{H}_{0}-\frac{i}{2}\sum_{k=1}^{L}\hat{L}^{k^{\dagger}}\hat{L}^{k}.$$

(2)

In following discussion, we omit subscript for the non-Hermitian Hamiltonian. Since the Lindbladian operator above, we are interested in exploring the properties of non-Hermitian many-body system which is described by the effective Hamiltonian as

$$\hat{H}=\hat{H}_{0}-\frac{i}{2}\Gamma\sum_{k=1}^{L}\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k}.$$

(3)

From above formula, which gives out the non-Hermitian Hamiltonian of a class of spin many-body system with jump operator of $\hat{L}^{k}=\sqrt{\Gamma}\sigma_{-}^{k}$, it can be found that dissipation imposes imaginary field along $z$-direction by looking at formula $\sigma_{+}^{k}\sigma_{-}^{k}=(\sigma_{z}^{k}+I)/2$. We look back at the knowledge about the general matrix. There are the left and right eigen-energy and eigenvectors with the non-Hermitian effective Hamiltonian, which are denoted by

	$\displaystyle\hat{H}|\phi_{R}^{i}\rangle=E_{i}|\phi_{R}^{i}\rangle,$		(4)
	$\displaystyle\langle\phi_{L}^{i}|\hat{H}=E_{i}\langle\phi_{L}^{i}|,$		(5)

with eigenvalues $E_{i}\in\mathbb{C}$. The left and right eigenvectors satisfy bi-orthonormal relation and completeness relation [19]

	$\displaystyle\langle\phi_{L}^{i}|\phi_{R}^{j}\rangle=\delta_{ij}\langle\phi_{L}^{j}|\phi_{R}^{j}\rangle,$		(6)
	$\displaystyle\sum_{j}\frac{|\phi_{R}^{j}\rangle\langle\phi_{L}^{j}|}{\langle\phi_{L}^{j}|\phi_{R}^{j}\rangle}=I.$		(7)

The properties of eigen-energy and eigenvectors of non-Hermitian are completely different from the Hermitian case, which endow non-Hermitian system with interesting new physics. When at least two eigen-energies coalesce with each other, exceptional point(EP) emerges[3, 4, 16, 19, 20, 21]. In order to obtain EP with our model in Eq(3), original Hamiltonian $H_{0}$ must not conserve spin polarization. If $\left[\hat{H}_{0},\sum_{k}\sigma_{z}^{k}\right]=0$ holds, the non-Hermitian term in Eq(3) can only give us trivial complex energy. Such that we choose quantum contact process(QCP) model in 1D as research object

$$\hat{H}_{0}=\Omega\sum_{k=1}^{L}\left(\hat{\sigma}_{x}^{k}\hat{\sigma}_{n}^{k+1}+\hat{\sigma}_{n}^{k}\hat{\sigma}_{x}^{k+1}\right),$$

(8)

with $\hat{\sigma}_{n}^{k}=\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k}$ and $L$ the number of spins in the system. QCP model has been widely studied in open many-body system, such as quantum epidemic dynamics [22, 23], experimental realization by Rydberg atom [24]. We arrive at non-Hermitian Hamiltonian as

$$\hat{H}=\Omega\sum_{k}^{L}\left(\hat{\sigma}_{x}^{k}\hat{\sigma}_{n}^{k+1}+\hat{\sigma}_{n}^{k}\hat{\sigma}_{x}^{k+1}\right)-\frac{i}{2}\Gamma\hat{\sum}_{k=1}^{L}\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k},$$

(9)

which determines non-Hermitian physics in our model. For convenience with following discussion, we set $\Omega=1$ to make the parameter of Hamiltonian dimensionless. In our work, we always use periodic boundary condition.

We firstly look at the complex energy spectrum of the non-Hermitian system for two spins case $L=2$. We obtain four eigen-energies as $E_{1}=0\quad(\rm trivial)$, $E_{2}=-i\Gamma/2$, $E_{3,4}=\frac{-i3\Gamma}{4}\pm\frac{1}{4}\sqrt{32\Omega^{2}-\Gamma^{2}}$, please see the [25] for the matrix and eigenvectors for $L=2$. It is easy to determine the EP as $\Gamma_{c}=\sqrt{32}\Omega$, where $E_{3}$ and $E_{4}$ coalesce with each other. For many-body case, we look at the energy spectrum changing with dissipation strength $\Gamma$ by exact diagonalization(ED). We used brute-force diagonalization up to $N=12$ as eig of matlab, and diagonalization based on Lanczos method for ground state up to $N=16$ as eigs of matlab. See Fig. 1 subplots(a-d), energy spectrum changes with $\Gamma$ from distributing on complex plane of $E$ to distributing on imaginary axis. Moreover, the real part of the energy varies with dissipation strength $\Gamma$ from finite values to zero at around $\Gamma\sim 14$, whose singularity looks like the order parameter of second-order phase transition [26]; while there is split with the imaginary part of energy. We also notice that $i\hat{H}_{eff}$ satisfies pseudo-Hermiticity [1] by transforming Hamiltonian and energy spectrum

	$\displaystyle\hat{H}_{eff}\rightarrow i\hat{H}_{eff}$
	$\displaystyle(\pm E_{re}+iE_{im})\rightarrow(\pm iE_{re}-E_{im})$

which results in that the eigenvalues appear in complex conjugate pairs at left of EP and the spectrum are entirely real at right of EP. These analysis show that $H_{eff}$ of NH QCP has anti-PT symmetry breaking [27, 28, 29] near EP. The singularity of the real part energy and anti-PT breaking indicates a transition in the energy spectrum of the non-Hermitian model (9), as a complex-imaginary transition.

In order to understand the complex-imaginary transition better, we are going to analyse the phase transition by studying the order parameter and susceptibility. But here we need to make preparation by defining the ground state of non-Hermitian many-body system. Ground states: we define the ”ground state” $|\psi_{R/L}\rangle$ as the state $\psi_{0}$ associated to the energy $E_{0}$ with minimum real part as did in literatures [14, 30]. The state with smallest real part can map to traditional ground state through continuation between the non-Hermitian case and Hermitian case, which can be observed within short time evolution in the absence of quantum jump[15, 14] due to its imaginary part. The state with largest imaginary part corresponds to the steady state that can be observed within long time evolution [21, 20]. By continuation we mean that the non-Hermitian Hamiltonian becomes the Hermitian case when the dissipation $\Gamma$ tends to zero and the chosen ”ground state” recover to the traditional one, namely the state corresponds to the smallest real energy. The former satisfies the requirement of continuation while the latter does not. When the system crosses to right side of the exceptional point, since energies become imaginary such that above definition of ground state become invalid, we track the continuity and analyticity of energy and observables to choose the ground state. Concretely, the ground state should keep the continuity of imaginary part of ground state from left side to right side of EP firstly. After reaching right side of EP, we use continuity and analyticity of observables to determine the ground state. In following part, we will reveal the phase transition induced by non-Hermiticity in our model in detail.

Generally, we use right ground eigenstate to compute the observables

$$O=\langle\psi_{R}|\hat{O}|\psi_{R}\rangle,$$

(10)

which can make sure that observables are always real. In fact, it can give the same results by using right eigenstate and using left-right eigenstate, see Fig.10 in appendix [25]. In order to compute the susceptibility of system, we impose a small probe magnetic field along the $z$-axis,

$$\hat{H}\rightarrow\hat{H}-\delta h/2\sum_{k}\hat{\sigma}_{z}^{k}=\hat{H}-\delta h\sum_{k}(\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k}-I/2).$$

(11)

The spin polarization $M^{\alpha}=\langle\psi_{R}|\sum_{k}\hat{\sigma}_{\alpha}^{k}|\psi_{R}\rangle$ for general direction $\alpha=x,y,z$. In our work, we choose spin polarization along $x$-axis $M^{x}$ as order parameter. The susceptibility can be given by the differential,

$$\chi=(M^{z}(\delta h)-M^{z}(0))/2\delta h.$$

(12)

The non-analytical behaviour of order parameter and susceptibility can help us to reveal the property of phase transition for non-Hermitian QCP model.

As shown in Fig. 2, order parameter $M^{x}$ vary with dissipation strength $\Gamma$ for different spin number up to $L=16$ in linear axis (a) and logarithmic axis (b). The critical points for different system size are determined when $M^{x}$ becomes zero. While in Fig. 2(b) the data almost collapse into one straight line in logarithmic scale, we extract the critical exponents $\beta$ for different spin number according to scaling function

$$M^{x}\sim-|\Gamma_{c}(L)-\Gamma|^{\beta}.$$

(13)

The critical points for different system size gradually converge to true critical point. In Table.1, we list critical points and critical exponents for system size from $L=4$ to $L=16$. And $\mathrm{normr}$ stands for norm of residuals, which quantifies the error for fitting our data to the scaling function in logarithmic axis. All critical exponents for different system sizes are nearly the same. After taking average of critical exponents of different $L$, the critical exponent of order parameter is estimated as $\bar{\beta}\approx 0.51$. In fact, we can look at theoretical results of special case $L=2$ to check numerical results. The order parameter can be written as

$$M^{x}\sim\sqrt{2\Gamma_{c}(2)}\sqrt{\Gamma_{c}(2)-\Gamma}\sim\sqrt{\Gamma_{c}(2)-\Gamma},$$

(14)

here $\Gamma_{c}(2)=\sqrt{32}\Omega$ mean the EP for $L=2$, see appendix [25] for details. This theoretical analysis confirm numerical results of critical exponent. From above analysis, we can confidently claim the true critical exponent for thermodynamic limit $\beta=1/2$ [31]. The results of order parameter $M^{x}$ confirm the analysis by energy spectrum in Fig. 1, namely there is continuous phase transition induced by non-Hermiticity in QCP model. To be added, we observe the non-analytical behaviour of order parameter for finite size system even $L=2$. While, for classical phase transition of Hermitian models, non-analytical behaviour of order parameter appears in thermodynamic limit. When system size is finite, order parameter become smooth near critical point. In reference [20], the authors noticed similar behaviours. The non-analytical behaviour of order parameter for finite system size in our work may be related to non-Hermiticity, which is very different from classical phase transition of Hermitian case.

We also study the magnetic susceptibility $\chi$ along $z$-direction, see Fig. 3. We show susceptibility vary with dissipation $\Gamma$ for different spin number up to $L=12$ in linear axis (a) and logarithmic axis (b). We can also observe that the susceptibility tends to infinite divergence for finite size system rather than finite peak for classical phase transition. As increasing system size, susceptibility strongly diverges near critical points in Fig. 3(a). Moreover these diverging data can collapse into one straight line in logarithmic scale in Fig. 3(b) except the data of small size $L=4$ with small error. The critical exponent of susceptibility is extracted according to scaling relation

$$\chi=\chi_{0}/(\Gamma-\Gamma_{c}(L))^{\gamma}.$$

(15)

In Table.2, we list critical points, critical exponents and amplitude of susceptibility for different system size. By comparing the critical points of Table.1 and Table.2, both results of $\Gamma_{c}(L))$ precisely agree with each other with negligible error. The critical exponents $\gamma$ are close to each other near $1.5$. After taking average of critical exponents $\gamma$, mean value of critical exponent is estimated as $\bar{\gamma}\approx 1.52$. Usually, EP is sensitive to disorder [32]. In appendix.C, we consider the non-Hermitian QCP model with random potentials and two spins case. The results indicate that the critical exponent changes even if under weak perturbation, the deviation from $\beta=0.5$ is slight within $5\%$ as shown in Fig. 17 of Appendix.C

In order to understanding the correlation and entanglement near phase transition, we study the correlation function of spin chain with open boundary

$$\Delta\sigma^{x}_{1}\sigma^{x}_{n}=\langle\hat{\sigma}_{x}^{1}\hat{\sigma}_{x}^{n}\rangle-\langle\hat{\sigma}_{x}^{1}\rangle\langle\hat{\sigma}_{x}^{n}\rangle.$$

(16)

we also look at Von Neumann entropy $S_{vn}$ of subsystem-$A$ by partitioning system as $L=L_{A}+L_{B}$

$$S_{vn}=-\Tr[\rho_{A}\ln\rho_{A}]=-\sum_{n}\lambda_{n}\ln\lambda_{n},$$

(17)

where reduced density matrix of $A$ is $\rho_{A}=\Tr_{B}[|\psi_{R}\rangle\langle\psi_{R}|]$ of ground state $|\psi_{R}\rangle$, and $\{\lambda_{n}\}$ are the eigenvalues of $\rho_{A}$. During calculation, we define the half-partition as

$$L_{A}=\left\{\begin{array}[]{ccc}L/2&,&\mod(L,2)=0\\ (L+1)/2&,&\mod(L,2)=1\\ \end{array}\right.$$

(18)

we also used the right eigenstate to make sure that correlation function and Von Neumann entropy are real numbers in the non-Hermitian model.

As seen in fig.4(a), the correlation functions vary with spin site $n$ for parameters $\Gamma=12,13$(red and blue lines) at left side and $\Gamma=14.5,15$(cyan and pink lines) right side of critical point. Correlations have quasi long-range decay at left side of critical point and nearly exponential decay at right side of critical point. The system at left side seem be quasi long-range ordering (see similar case [20]). Meanwhile, system at right side has short-range ordering. The correlation at left and right side of critical point decay in completely different form regardless of boundary condition, as fig.13 in appendix left panel for periodic boundary and right for open boundary .

As seen in Fig.4(b), entanglement entropy decays with $\Gamma$ for different $L$. Non-analytical peaks indicate a NH phase transition for the system from a higher to a lower entangled state. For ground state of NH-QCP model, entanglement characterize the continuous phase transition by non-analytical peak. At same time, we consider two cases to study entanglement entropy, first one: $L$ changes from $4$ to $16$ for half-partition $L_{A}=L/2$ and second one: $L_{A}$ changes from $1$ to $L/2$ for fixed $L=16$. In Fig.5, we also found the results that entanglement entropy $S_{vn}$ almost (a) do not change with $L$ for half- partition and (b) do not vary with $L_{A}$ for $L_{A}$ ranging from $1$ to $L/2$ when $\Gamma\in[12.0,13.0,14,15.0]$ is near $\Gamma_{c}$. These results reveal that entanglement entropy in NH-QCP model satisfies area law since the surface of spin chain is two spins as Hermitian model in 1D.

By summarizing above results, we have obtained the critical exponents $\beta$, $\gamma$ for the phase transition in non-Hermitian QCP model. As far as we know, $\beta$ in our model is same as mean field results; While $\gamma$ in our work are different from the known critical exponents. Since the phase transition is related with the coalescence of energy spectrum at EP, non-Hermiticity can endow non-Hermitian many-body system with new critical behaviour, such as shown in Fig. 2. 3. The analysis about critical exponents of $\beta$, $\gamma$ reveal that the phase transition in non-Hermitian QCP is independent of system size. But the critical points serve as function of system size $L$, and correlation length must be constrained by $L$. The critical behaviour is governed by the long-wavelength fluctuations such that the critical exponents are universal. For finite size, singularity of phase transition in this work accompany with finite correlation length. This peculiar behaviour should be derived from the singular behaviour of EP for non-Hermitian system. The phase transition for non-Hermitian many-body system may be beyond the classification for Hermitian many-body system.

We compare the results of our model with its Hermitian counterpart by generalizing parameter to complex value $i\Gamma\rightarrow\Gamma_{\rm re}+i\Gamma_{\rm im}$ in Eq.(9). The Hamiltonian becomes

$$\hat{H}=\Omega\sum_{k}\left(\hat{\sigma}_{x}^{k}\hat{\sigma}_{n}^{k+1}+\hat{\sigma}_{n}^{k}\hat{\sigma}_{x}^{k+1}\right)-\frac{1}{2}(\Gamma_{\rm re}+i\Gamma_{\rm im})\sum_{k=1}^{L}\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k},$$

(19)

where we define $\Gamma=\Gamma_{\rm re}+i\Gamma_{\rm im}$. When $\Gamma$ taking real part $\Gamma_{\rm re}$, it corresponds to Hermitian case; While $\Gamma$ taking imaginary part $i\Gamma_{\rm im}$, it corresponds to non-Hermitian case in this work. In Fig. 6, we show the susceptibility of non-Hermitian case (a) and Hermitian counterpart(b). By numerical fitting, the critical exponents $\gamma$ for both case are respectively near $1.5,0.57$. From the difference of critical exponent $\gamma$, it indicates that the continuous phase transition in non-Hermitian QCP model is different from its Hermitian counterpart. By looking at the number of up-spin(c)and order parameter(d) for Hermitian case, these observables discontinuously jump at critical point $\Gamma_{\rm re}^{c}\sim-4$, which imply first-order phase transition with Hermitian QCP model.

We include complementary results of Hermitian model in Figs.7.8. As Fig.7 of numerical results shows, the lowest five energy levels vary with $\Gamma_{\rm re}$. There is level crossing point at $\Gamma_{\rm re}=-4\Omega$ for Hermitian QCP. As shown in Fig.8, the number of up-spin and order parameter which are rescaled with $L$ vary with $\Gamma_{\rm re}$ for different system sizes, indicating that the jump heights are proportional to system size $L$. There are sudden jumps with both observables near critical point. All results for different system size collapse into same curve.

Symmetry analysis: We transform the Hamiltonian of Hermitian case by using $\sigma_{n}^{k}=(\sigma_{z}^{k}+1)/2$

		$\displaystyle\hat{H}$	$\displaystyle=\frac{\Omega}{2}\sum_{k}\left(\hat{\sigma}_{x}^{k}\hat{\sigma}_{z}^{k+1}+\hat{\sigma}_{z}^{k}\hat{\sigma}_{x}^{k+1}\right)+\Omega\sum_{k}\hat{\sigma}_{x}^{k}$		(20)
		$\displaystyle-$	$\displaystyle\frac{1}{4}\Gamma_{\rm re}\sum_{k=1}^{L}\hat{\sigma}_{z}^{k}+const.$

We notice that in the interaction terms the $x$ and $z$ directions play a symmetric role. On the other hand, the external fields (second and third terms) break this symmetry. The symmetric point is $\Gamma_{\rm re}=-4\Omega$, which agrees with critical point in Figs.7.8 by numerics. We can conclude that order parameter $M^{x}$ of ground state has sudden jump at two side of critical point hence first-order phase transition with Hermitian counterpart.

In fact, by summarizing the results of critical exponent, we think that phase transition in non-Hermitian QCP model belongs to new universality class from previous known’s. Non-Hermitian many-body system break the Hermitian conjugation symmetry. This type of phase transition is induced by the non-Hermiticity. Moreover, exceptional point is continuous and non-analytical behaviour mathematically for finite size system even $L=2$. At EP, non-Hermiticity alters the property of ground state hence the non-analytical behaviour of order parameter and singularity of susceptibility. It is uniqueness of non-Hermitian system comparing with regular phase transition of Hermitian counterpart.

We have studied the phase transition in non-Hermitian QCP model. We have determined the critical exponents $\beta$, $\gamma$ which may indicate that the phase transition in non-Hermitian QCP model belongs to new universality class. Our results can help people to understand the property of non-Hermitian many-body system. Moreover, our results can stimulate more study about the phase transition of non-Hermitian many-body system. As for outlook, the dynamics of non-Hermitian system will be interesting topic since non-unitary evolution is completely different from the unitary evolution of Hermitian model, like [33, 34]. Renormalization group can be used to reveal the difference and relation of phase transition between the Hermitian model and non-Hermitian model [35]. It may be feasible to utilize Rydberg atom[24, 36] to simulate the non-Hermitian QCP model in our work to observe the continuous phase transition experimentally even for simple case of two spins.

We gratefully thank Profs.Rosario Fazio and Stefano Chesi for their inspiring discussions and suggestions. W.B.H. acknowledges support from NSAF U1930402. J.J. acknowledges support from the National Natural Science Foundation of China (NSFC) via Grant No. 11975064. F.I. acknowledges the financial support of the Brazilian funding agencies National Council for Scientific and Technological Development-CNPq(Grant No. 308205/2019 − 7) and FAPERJ (Grant No. E-26/211.318/2019 and E-26/201.365/2022). H.Q. Lin acknowledges financial support from NSAF U1930402 and NSFC 12088101 and computational resources from the Beijing Computational Science Research Center.

The data that support the findings of this study are available upon reasonable request from the authors.