††thanks: These two authors contributed equally to the work††thanks: These two
authors contributed equally to the work

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

Wen-Bin He  Beijing Computational Science Research Center, Beijing 100193,
China  Quantum Systems Unit, Okinawa Institute of Science and Technology
Graduate University, 904-0495 Okinawa, Japan  Jiasen Jin  School of Physics,
Dalian University of Technology, Dalian 116024, China  Fernando Iemini
<EMAIL_ADDRESS>Instituto de Física, Universidade Federal Fluminense,
24210-346 Niterói, Brazil  The Abdus Salam International Center for
Theoretical Physics, Strada Costiera 11, 34151 Trieste, Italy.  Hai-Qing Lin
<EMAIL_ADDRESS>Beijing Computational Science Research Center, Beijing
100193, China  Department of Physics, Beijing Normal University, Beijing
100875, China

###### Abstract

Non-Hermitian quantum system recently have attracted a lots of attentions
theoretically and experimentally. However, the results based on the single-
particle picture may not apply to understand the property of non-Hermitian
many-body system. How the property of quantum many-body system especially the
phase transition will be affected by the non-hermiticity remains unclear. Here
we study non-Hermitian quantum contact process (QCP) model, whose effective
Hamiltonian is derived from Lindbladian master equation. We show that there is
a continuous phase transition induced by the non-hermiticity in QCP. We also
determine the critical exponents $\beta$ of order parameter, $\gamma$ of
susceptibility and study the correlation and entanglement near phase
transition. We observe that the order parameter and susceptibility display
infinitely singularity even for finite size system, since non-hermiticity
endow many-body system with different singular behaviour from classical phase
transition. Moreover our results show that the phase transition have no
counterpart in Hermitian case and belongs to completely different universality
class.

Introduction.  Non-Hermitian physics currently has attracted considerable
attentions Ueda_review ; Bergholtz , which greatly extend the knowledge about
quantum system. For example, the complex energy coalesce at exceptional
point(EP)Bergholtz ; Hatano ; wang18 ; The wave function localize at the edge
of quantum system, which is now named as non-Hermitian skin effectHatano ;
wang18 ; Gong20 ; PXue18 . 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 Clara_nc ; Clara_l . 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
wang18 ; ZhongWang ; kazuki . Non-hermiticity can also induce phase transition
without gapNorifumi  and non-Hermitian many-body localization Hamazaki19 .
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$, correltion, and entanglement
entropy. The results of this work show that there is continuous phase phase
with non-hermitian QCP model. We extract the critical exponents, $\beta$ of
order parameter, $\gamma$ of susceptibility. And we reveal the transition from
quasi long-rang order to short-range order and entanglement effect 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.

Non-Hermitian QCP model. For the open quantum many-body systems, its dynamics
usually can be determined by Lindbladian master equation breuer  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 induces local spin decays. The
Hamiltonian for non-Hermitian system is associated with 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 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}$ , dissipation impose 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}$. And the left and right eigenvectors
satisfy bi-orthonormal relation and completeness relation Heiss

$\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) emergesHatano ; wang18 ; breuer ; Heiss ;
Lee ; Tony . 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, like quantum epidemic dynamics Lesanovsky ; Wintermantel ,
experimental realization by Rydberg atom Gutierrez . 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 determine 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.

Figure 1: The subplots (a)(b)(c)(d), distribution pattern of complex energy
spectrum for four different $\Gamma$; The lowest and highest seven real part
(e) and their imaginary part (f) of energy spectrum $E$ vary with $\Gamma$.
$\Omega=1,L=10$.

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 SMsm  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), see
Fig. 1. From subplots(a-d), energy spectrum change with $\Gamma$ from
distributing on complex plane of $E$ to distributing on imaginary axis.
Moreover, the real part of energy vary with dissipation strength $\Gamma$ from
finite values to zero at around $\Gamma\sim 14$; while there is split with the
imaginary part of energy. These results indicate a transition of energy
spectrum of non-Hermitian model (9), as complex-imaginary transition, like
second-order phase transition Sachdev . We also notice that $i\hat{H}_{eff}$
satisfies pseudo-hermiticity Ueda_review  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 Ramy ; ozdemir ; YYang .

In order to understand the complex-imaginary transition better, we are going
to analyze 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 Norifumi ; weitao . With
this definition, by considering Hermitian limit ($\Gamma\rightarrow 0)$, we
recover the standard ground state definition. 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. Concretly,
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.

Figure 2: (a)The spin polarization along $x$-axis direction $M^{x}=\left\langle\sum_{k}\sigma^{k}_{x}\right\rangle$ of ground state vary with $\Gamma$ for different system size $L$; (b) the absolute value of $M^{x}$ near critical point V.S dissipation strength in logarithmic scale. The scaling exponent of $M^{x}\sim-|\Gamma_{c}(L)-\Gamma|^{\beta}$ is fitted out for different system size, see the Table.(1). $\Omega=1$.  Table 1: The critical points,critical exponents of order parameter of spin polarization along x-direction $M^{x}$ for different system size $L$. $\mathrm{normr}$ stands for norm of residuals. $L$ |  4 |  5 |  6 |  7 |  8 |  9 |  10
---|---|---|---|---|---|---|---
$\Gamma{c}(L)$ |  13.4396 |  13.6116 |  13.6969 |  13.7424 |  13.7720 |  13.7898 |  13.8016
$\beta$ |  0.5102 |  0.5209 |  0.5303 |  0.5117 |  0.5261 |  0.5251 |  0.5026
$\mathrm{normr}$ |  0.0336 |  0.0582 |  0.0953 |  0.0108 |  0.0319 |  0.0383 |  0.0122
$L$ |  11 |  12 |  13 |  14 |  15 |  16 |
$\Gamma{c}(L)$ |  13.8106 |  13.8169 |  13.8223 |  13.8259 |  13.8289 |  13.8312 |
$\beta$ |  0.5069 |  0.4968 |  0.5178 |  0.5120 |  0.5201 |  0.5122 |
$\mathrm{normr}$ |  0.0111 |  0.0687 |  0.0063 |  0.0129 |  0.0504 |  0.0181 |

The observables and critical exponents. 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. 1 in SM sm . 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 almostly
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 quantify the error for fitting our data to the scaling
function in logarithmic axis. All critical exponents for different system size
nearly are 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 SMsm  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$ BoWei . 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 Lee , 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.

Table 2: The critical points, critical exponents and amplitude of susceptibility for different system size $L$. $L$ |  4 |  5 |  6 |  7 |  8
---|---|---|---|---|---
$\Gamma{c}(L)$ |  13.4388 |  13.6092 |  13.6949 |  13.7418 |  13.7704
$\gamma$ |  1.5198 |  1.5829 |  1.4975 |  1.5049 |  1.4986
$\chi_{0}*10^{-4}$ |  0.2972 |  0.2614 |  0.3377 |  0.3384 |  0.3573
$\mathrm{normr}$ |  0.0009 |  0.0061 |  0.0001 |  0.0001 |  0.0003
$L$ |  9 |  10 |  11 |  12 |
$\Gamma{c}(L)$ |  13.7887 |  13.8010 |  13.8105 |  13.8168 |
$\gamma$ |  1.5273 |  1.4969 |  1.4999 |  1.5731 |
$\chi_{0}*10^{-4}$ |  0.3337 |  0.3907 |  0.3809 |  0.2833 |
$\mathrm{normr}$ |  0.0008 |  0.0022 |  0.0002 |  0.0110 |
Figure 3:  The susceptibility $\chi$ vary with $\Gamma$ in linear scale(a) and
vary with $\Gamma-\Gamma_{c}(L)$ in logarithmic scale(b) for different size
$L=[4,5,6,7,8,9,10,11,12]$ . By linear fitting in logarithmic scaling of
formula $\chi=\chi_{0}/(\Gamma-\Gamma_{c}(L))^{\gamma}$, the critical points,
the scaling exponents $\gamma$ and amplitude $\chi_{0}$ in the table.2.

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 susceptibility tend to infinite divergence for finite size system
rather than finite peak for classical phase transition. As increasing system
size, susceptibility strongly diverge 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. And 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$. These results give us strong clue that the true critical exponent in
thermodynamic limit $\gamma$ of susceptibility should be $1.5$.

Figure 4:  The correlation function $\Delta\sigma^{x}_{1}\sigma^{x}_{n}$ (a)
vary with spin site $n$ for four different $\Gamma$ with system size $L=16$.
Von Neumann entropy $S_{vn}$(b) vary with $\Gamma$ for half-partition of
different system size $L$.

Correlation and entanglement. 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)

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

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

where reduced density matrix of $A$ is
$\rho_{A}=\operatorname{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 consider half-partition $L_{A}=L/2$.

As seen in fig.4(a), the correlation function 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.
Correlation 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 Lee ). Meanwhile,
system at right side has short-range order. 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.
In SMsm  figs.7.8, we also include the results that $S_{vn}$ almost do not
change with $L$ for half- partition and $L_{A}$-partition for $L_{A}$ ranging
from $1$ to $L/2$ when $\Gamma\in[12.0,13.0,14,15.0]$ is away from
$\Gamma_{c}$. These results reveal that entanglement entropy in NH QCP model
satisfy 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 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.

_Comparison with Hermitian counterpart._ 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},$
(18)

and 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. 5, 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}\approx-4$, which imply first-order phase transition with
Hermitian QCP model.

We include complementary results of Hermitian model in SMsm , which contain
the energy, the number of up-spin, and order parameter $M^{x}$. For
counterpart of Hermitian QCP, there is energy level crossing with ground
state, meanwhile the number of up-spin and order parameter $M^{x}$ have
suddenly jump near the critical point and collapse onto the same curve for
different system size. The interested reader can find more details in SMsm .
These results indicate that the phase transition for Hermitian model is first-
order. 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 case.

Figure 5: subplot(a): By setting real part of $\Gamma$ as zeros, the
susceptibility vary with $\Gamma_{im}$ for ground eigenstates of non-Hermitian
QCP model. Subplot(b):By setting imaginary part of $\Gamma$ as zeros, the
susceptibility of Hermitian QCP model vary with $\Gamma_{re}$ for ground
eigenstates. The critical exponent for Hermitian model is $\gamma\sim 0.5699$.
Subplots(c)(d). The $n_{up}$ number of up-spin and order parameter $M^{x}$
vary with $\Gamma_{re}$. $L=12$, $\Omega=1$.

Conclusion. 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 Haga ; Orito . Renormalization group can be used to
reveal the difference and relation of phase transition between the Hermitian
model and non-Hermitian model Kawabata1 . It may be feasible to utilize
Rydberg atomGutierrez ; Browaeys  to simulate the non-Hermitian QCP model in
our work to observe the continuous phase transition experimentally even for
simple case of two spins.

Acknowledgements. We gratefully thank Prof.Rosario Fazio and Stefano Chesi for
their inspiring discussions and suggestions. W.B.H acknowledges support from
Grant No. U1930402 of the National Science Association Funds of NSFC. 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$). H.Q.L. acknowledges financial support from National
Science Association Funds U1930402 and NSFC 11734002, as well as computational
resources from the Beijing Computational Science Research Center.

## References

  * (1) Yuto Ashida, Zongping Gong, and Masahito Ueda, Non-Hermitian physics, Advances in Physics 69 249 (2020).
  * (2) Emil J. Bergholtz, Jan Carl Budich, and Flore K. Kunst, Rev. Mod. Phys. 93, 015005(2021).
  * (3) N. Hatano and D. R. Nelson, Phys. Rev. Lett. 77, 570 (1996).
  * (4) Shunyu Yao and Zhong Wang, Phys. Rev. Lett. 121, 086803(2018).
  * (5) Linhu Li, Ching Hua Lee, Sen Mu and Jiangbin Gong, Nat.Comm.11,5491 (2020).
  * (6) L. Xiao, T. Deng, K. Wang, G. Zhu, Z. Wang, W. Yi, and P. Xue, Non-Hermitian bulk-boundary correspondence in quantum dynamics, Nat. Phys. 16, 761 (2020).
  * (7) Clara C. Wanjura, Matteo Brunelli,and Andreas Nunnenkamp Topological framework for directional amplification in driven-dissipative cavity arrays, Nat. Comm. 11, 3149 (2020).
  * (8) Clara C. Wanjura, Matteo Brunelli, and Andreas Nunnenkamp Topological framework for directional amplification in driven-dissipative cavity arrays, Phys.Rev.Lett. 127,213601(2021).
  * (9) Shunyu Yao, Fei Song, and Zhong Wang, Non-Hermitian Chern Bands, Phys. Rev. Lett.121,136802(2018).
  * (10) Kazuki Yokomizo and Shuichi Murakami, Non-Bloch Band Theory of Non-Hermitian Systems, Phys. Rev. Lett. 123, 066404(2019).
  * (11) Norifumi Matsumoto, Kohei Kawabata, Yuto Ashida, Shunsuke Furukawa, and Masahito Ueda, Continuous Phase Transition without Gap Closing in Non-Hermitian Quantum Many-Body Systems, Phys. Rev. Lett. 125, 260601(2020).
  * (12) Ryusuke Hamazaki, Kohei Kawabata, and Masahito Ueda,Phys. Rev. Lett. 123, 090603(2019).
  * (13) H. P. Breuer, and F. Petruccione, Theory of Open Quantum Systems, Oxford University Press, Oxford, (2002).
  * (14) W D Heiss The physics of exceptional points, J. Phys. A: Math. Theor. 45, 444016(2012).
  * (15) Tony E. Lee and Ching-Kit Chan, Heralded Magnetism in Non-Hermitian Atomic Systems, Phys.Rev.X.4, 041001(2014).
  * (16) Tony E. Lee, Florentin Reiter, and Nimrod Moiseyev,Entanglement and Spin Squeezing in Non-Hermitian Phase Transitions, Phys. Rev. Lett. 113, 250401(2014).
  * (17) Federico Carollo, Edward Gillman, Hendrik Weimer, and Igor Lesanovsky,Critical Behavior of the Quantum Contact Process in One Dimension, Phys. Rev. Lett. 123, 100604 (2019).
  * (18) T. M. Wintermantel, M. Buchhold, S. Shevate, M. Morgado, Y. Wang, G. Lochead, S. Diehl, and S. Whitlock,Epidemic growth and Griffiths effects on an emergent network of excited atoms Nat Commun 12 103 (2021).
  * (19) Ricardo Gutiérrez, Cristiano Simonelli, Matteo Archimi, Francesco Castellucci, Ennio Arimondo, Donatella Ciampini, Matteo Marcuzzi, Igor Lesanovsky, and Oliver Morsch, Phys. Rev. A 96, 041602(R) (2017).
  * (20) Supplementary material for more details.
  * (21) S. Sachdev, Quantum Phase Transitions (Cambridge University Press, Cambridge, England, 2001).
  * (22) Ramy El-Ganainy, Konstantinos G. Makris, Mercedeh Khajavikhan, Ziad H. Musslimani, Stefan Rotter and Demetrios N. Christodoulides, Non-Hermitian physics and PT symmetry, Nat. Phys. 14, 11–19 (2018).
  * (23) S.K. Özdemir, S. Rotter, F. Nori and L. Yang , Parity–time symmetry and exceptional points in photonics, Nature Materials. 18, 783–798 (2019).
  * (24) Y. Yang, Yi-Pu Wang, J. W. Rao, Y. S. Gui, B. M. Yao, W. Lu, and C.-M. Hu, Unconventional Singularity in Anti-Parity-Time Symmetric Cavity Magnonics, Phys. Rev. Lett. 125, 147202(2020).
  * (25) Weitao Chen, Liangtao Peng, Hantao Lu, Xiancong Lu, Phys. Rev. B 105, 075126(2022).
  * (26) Bo-Bo Wei and Liang Jin, Universal Critical Behaviours in Non-Hermitian Phase Transitions, Scientific Reports. 7,7165 (2017) .
  * (27) Taiki Haga, Masaya Nakagawa, Ryusuke Hamazaki, and Masahito Ueda,Liouvillian Skin Effect: Slowing Down of Relaxation Processes without Gap Closing, Phys. Rev. Lett. 127, 070402(2021).
  * (28) Takahiro Orito, Ken-Ichiro Imura,Unusual diffusion and entanglement dynamics in non-Hermitian disordered many-body systems,Phys. Rev. B 105, 024303(2022) .
  * (29) Kohei Kawabata1, and Shinsei Ryu,Nonunitary Scaling Theory of Non-Hermitian Localization, Phys.Rev.Lett.126,166801(2021).
  * (30) A. Browaeys and T. Lahaye, Many-body physics with individually controlled Rydberg atoms, Nat. Phys. 16, 132 (2020)

##  I  Supplementary Material: The continuous phase transition induced by the
non-Hermiticity in quantum contact process model.

Wen-Bin He Jiasen Jin Fernando Iemini

Hai-Qing Lin

In this supplementary material(SM), we include more results which can help
understand the main text of our work. The SM contains the results of numerical
results for general $L$, special case $L=2$, Hermitian case.

##  II The definitions and notations

In this part, we firstly introduce the definition used in main text and show
more numerical results.

Observables.- Generally, we use right eigenstate of ground state to compute
the observables:

$O=\langle\psi_{R}|\hat{O}|\psi_{R}\rangle.$ (19)

such as, number of spin-up

$n_{\rm
up}=\langle\psi_{R}|\sum_{k}\hat{\sigma}_{+}^{k}\hat{\sigma}_{-}^{k}|\psi_{R}\rangle.$
(20)

Susceptibility.-  In order to compute the susceptibility, we impose a small
probe magnetic field along the $z$-axis,

$\hat{H}\rightarrow\hat{H}-\delta h\sum_{k}\hat{\sigma}_{z}^{k}.$ (21)

The susceptibility can be given as,

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

with magnetization
$M^{\alpha}=\langle\psi_{R}|\sum_{k}\hat{\sigma}_{\alpha}^{k}/2|\psi_{R}\rangle$,
for $\alpha=x,y,z$. In this way, the observables and susceptibility are always
real which is conformal to experimental detection.

We can also use the left and right eigenstates of ground state to define the
observable

$O_{LR}=\langle\psi_{L}|\hat{O}|\psi_{R}\rangle,$ (23)

and its susceptibility is denoted as $\chi_{LR}$ respectively. We think that
$\chi$ and $\chi_{LR}$ both have the same function to characterize the
property of phase transitions, see fig.6. Near critical point,$\chi$ is very
to close to the real part of $\chi_{LR}$. In figs.7.8, we show the order
parameter for $L=16$ and susceptibility up to $L=12$, which strongly indicate
a continuous phase transition.

Energy gap.  We define the energy gap as following,

$\Delta=\min_{n\neq 0}|E_{n}-E_{0}|.$ (24)

In fig.9, the energy gap for no-Hermitian Hamiltonian vary with dissipation
strength for $L=10$, which indicates the gap closing at critical point.

In fig.10, we show the extrapolation of the critical points $\Gamma_{c}$ and
the system size $L$. According to numerical results, the critical points of
thermaldynamic limit $L\rightarrow\infty$ is $\Gamma_{c}\sim 13.845$.

Figure 6: subfig(a):The susceptibility vary with $\Gamma$ for right ground
eigenstates of non-Hermitian QCP model. subfig(b):The susceptibility vary with
$\Gamma$ for right and left ground eigenstates of non-Hermitian QCP model.
$L=10$, $\Omega=1$. Figure 7: The spin polarization along $x$-axis direction
$M^{x}=\left\langle\sum_{k}\sigma^{k}_{x}\right\rangle$ vary with $\Gamma$ for
ground state of different system size $L=16$, $\Omega=1$. the scaling exponent
of $\beta=0.5122$, $\Gamma_{c}(L)=13.8312$ Figure 8:  The susceptibility vary
with $\Gamma$ for ground eigenstates for different $L$, $\Omega=1$. The black
dash lines indicate the critical points $\Gamma_{c}(L)$ of each spin number,
which are $13.4388,13.6949,13.770,13.801,13.8168$. We assume susceptibility
satisfy scaling law $\chi\sim 1/(\Gamma-\Gamma_{c})^{\gamma}$. By using the
right data of $\Gamma_{c}(L)$ which are larger than the critical point of each
cluster size, the critical exponents of each cluster can be extracted.  Figure
9: The energy gap $\Delta$ vary with $\Gamma$ for ground state $L=10$,
$\Omega=1$.

##  III analytical results of $L=2$

We look at the particular and simpler case of $L=2$ spins, which can provide
theoretical information about phase transition. It is easy to explicitly write
matrix of the effective Hamiltonian

$H_{\rm eff}=\left(\begin{array}[]{cccc}-i\Gamma&\Omega&\Omega&0\\\
\Omega&-i/2\Gamma&0&0\\\ \Omega&0&-i/2\Gamma&0\\\ 0&0&0&0\\\
\end{array}\right).$ (25)

The eigenvalues are given by

$\displaystyle E_{1}=0\quad(\rm trivial),$ $\displaystyle E_{2}=-i\Gamma/2,$
$\displaystyle
E_{3,4}=\frac{-i3\Gamma}{4}\pm\frac{1}{4}\sqrt{32\Omega^{2}-\Gamma^{2}}.$

Here, we only give out right eigenstates $|\psi_{R}^{i}\rangle$ since we use
right eigenstates to compute the observables

$\displaystyle|\psi_{R}^{1}\rangle$ $\displaystyle=$
$\displaystyle\left(0,0,0,1\right)^{T},$ $\displaystyle|\psi_{R}^{2}\rangle$
$\displaystyle=$ $\displaystyle\left(0,1,-1,0\right)^{T},$ (26)
$\displaystyle|\psi_{R}^{3,4}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2\Omega^{2}+\Gamma^{2}+|\lambda|^{2}+2\Gamma
Im(\lambda)}}\left(\sqrt{2}\Omega,(\lambda+i\Gamma)/\sqrt{2},(\lambda+i\Gamma)/\sqrt{2},0\right)^{T}.$

where, for last two wave functions, we remark that eigenvalue $\lambda$ are
$E_{3,4}$. The exceptional point(EP) is located at

$\Gamma_{c}=\sqrt{32}\Omega.$ (27)

When the dissipation $\Gamma$ exceed the $\Gamma_{c}$, all non-trivial
eigenvalues $E$ become the purely imaginary number. Meanwhile, the ground
state is $|\psi_{R}^{4}\rangle$ since $E_{4}$ has minimum real part.

The order parameter:
$M^{x}=\langle\psi_{R}^{4}|\sigma^{x}|\psi_{R}^{4}\rangle$, where
$\sigma^{x}=\sum_{k}\sigma^{k}_{x}$ can be written as matrix form

$\sigma^{x}=\left(\begin{array}[]{cccc}0&1&1&0\\\ 1&0&0&1\\\ 1&0&0&1\\\
0&1&1&0\\\ \end{array}\right).$ (28)

We obtain the results of $M^{x}$

$M^{x}=\frac{-\Omega\sqrt{32\Omega^{2}-\Gamma^{2}}}{2\Omega^{2}+\Gamma^{2}+|\lambda|^{2}+2\Gamma
Im(\lambda)}.$ (29)

When dissipation tends to exceptional point $|\Gamma-\Gamma_{c}|\rightarrow
0$, the order parameter can be written as
$M^{x}\sim\sqrt{2\Gamma_{c}}\sqrt{\Gamma_{c}-\Gamma}\sim\sqrt{\Gamma_{c}-\Gamma}$.
Such that the critical exponent $\beta=1/2$

Figure 10:  The extrapolation of critical points $\Gamma_{c}(L)$ of order
parameter $M^{x}$ for system size up to $L=16$. We extract the true critical
point $\Gamma_{c}\sim 13.845$.

##  IV  correlation and entanglement

Correlation function we used in this work is defined as

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

The correlation at left and right side of critical point decay in completely
different form(fig.6 left for periodic boundary, fig.6 right for open boundary
). By looking at the results of open boundary, at left side of critical point,
correlation seem be quasi long-range ordering. Meanwhile, correlation at right
side of EP seem be short-range.

Meanwhile we are also interested in the entanglement of phase transition. We
divide system by partition as $L=L_{A}+L_{B}$ , and study the entanglement
entropy of subsystem-$A$ as

$S_{vn}=-\operatorname{Tr}[\rho_{A}\ln\rho_{A}]=-\sum_{n}\lambda_{n}\ln\lambda_{n}$
(31)

where reduced density matrix of $A$ is
$\rho_{A}=\operatorname{Tr}_{B}[|\psi_{0}\rangle\langle\psi_{0}|]$ of ground
state $|\psi_{0}\rangle$. And $\\{\lambda_{n}\\}$ are the eigenvalues of
$\rho_{A}$. We define half-partition

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

For ground state of NH-QCP model, entanglement characterize the continuous
phase transition by non-analytical jump, see fig.4(b) in main text. 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$. Entanglement entropy
$S_{vn}$ almost do not change with system $L$ and $L_{A}$ for both cases, when
$\Gamma$ is away from $\Gamma_{c}$, see fig.13.14. These results indicate that
entanglement entropy satisfy area law as the surface of NH-QCP spin chain is
two spins.

Figure 11: Correlation function for periodic(left) and open(right) boundary.
$L=12$. The correlation length for $\Gamma=14.5$,$\Gamma=15.0$ are
$0.5879$,$0.5304$. Figure 12: Left panel: correlation function for open
boundary, $L=16$. For $\Gamma$ range from near $14$ to $16$ of right side of
critical point. Right panel: correlation length vary with $\Gamma$. $\nu$ is
small as $0.1586$(the data close to critical point give smaller $\nu$). .
Figure 13: Entanglement entropy $S_{vn}$ for half-partition vary with system
size $L$ for different $\Gamma=12.0,13.0,14,15.0$.  Figure 14: Entanglement
entropy $S_{vn}$ vary with $L_{A}$-partition $L=16$.
$\Gamma=12.0,13.0,14,15.0$.  Figure 15: The lowest five energy levels vary
with real part of parameter $\Gamma_{re}$ for hermitian case. $L=16$,
$\Omega=1$. Figure 16: The results for Hermitian model. subfig(a):The rescaled
$n_{up}$ vary with $\Gamma_{re}$ for ground eigenstates of Hermitian QCP model
for different system size $L$. subfig(b): rescaled $M^{x}$ vary with
$\Gamma_{re}$ for different system size $L$. $\Omega=1$.

##  V Hermitian cases

By combining the non-Hermitian QCP model and its Hermitian counterpart, we
obtain below Hamiltonian

$\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}.$
(33)

The complex-field $\Gamma=\Gamma_{\rm re}+i\Gamma_{\rm im}$ consists of real
part $\Gamma_{\rm re}$ and imaginary part $\Gamma_{\rm im}$, which have the
physical meaning as magnetic field along z-axis and dissipation respectively.
In this part, we particularly consider Hermitian case, namely $\Gamma_{\rm
im}=0$ for special case $L=2$.

Analysis of special case $L=2$.-  It is easy to explictly write the
Hamiltonian for this particular case Eq.(33)

$H=\left(\begin{array}[]{cccc}-\Gamma&\Omega&\Omega&0\\\
\Omega&-1/2\Gamma&0&0\\\ \Omega&0&-1/2\Gamma&0\\\ 0&0&0&0\\\
\end{array}\right).$ (34)

The eigenvalues are given by

$\displaystyle E_{1}=0\quad(\rm trivial),$ $\displaystyle E_{2}=-\Gamma/2,$
$\displaystyle
E_{3,4}=\frac{-3\Gamma}{4}\pm\frac{1}{4}\sqrt{32\Omega^{2}+\Gamma^{2}}.$

If there is an exceptional point,
$32\Omega^{2}+\Gamma^{2}=32\Omega^{2}+\Gamma_{re}^{2}-\Gamma_{im}^{2}+i2\Gamma_{re}\Gamma_{im}=0$.
We can solve it

$32\Omega^{2}+\Gamma_{re}^{2}-\Gamma_{im}^{2}=0,\quad\&\quad\Gamma_{re}\Gamma_{im}=0.$
(35)

Only for $\Gamma_{re}=0$ there is a reasonable solution, but recovers the
results of Eq.(9). Through the above analysis, for general non-zero
$\Gamma_{\rm re}$ and $\Gamma_{\rm im}$ there are no exceptional points of
energy coalescence.

If $\Gamma$ is real value, the eigenvalues

$\displaystyle E_{1}=0\quad(\rm trivial),$ $\displaystyle
E_{2}=-\Gamma_{re}/2,$ $\displaystyle
E_{3,4}=\frac{-3\Gamma_{re}}{4}\pm\frac{1}{4}\sqrt{32\Omega^{2}+\Gamma_{re}^{2}}.$

eigen-states are given as

$\displaystyle|\psi_{R}^{1}\rangle$ $\displaystyle=$
$\displaystyle\left(0,0,0,1\right)^{T},$ $\displaystyle|\psi_{R}^{2}\rangle$
$\displaystyle=$ $\displaystyle\left(0,1,-1,0\right)^{T},$ (36)
$\displaystyle|\psi_{R}^{3,4}\rangle$ $\displaystyle=$
$\displaystyle\frac{1}{\sqrt{2\Omega^{2}+(E+\Gamma_{\rm
re})^{2}}}\left(\sqrt{2}\Omega,(E+\Gamma_{\rm re})/\sqrt{2},(E+\Gamma_{\rm
re})/\sqrt{2},0\right)^{T}.$

Where $E$ takes $E_{3,4}$ for real $\Gamma$. If $\Gamma$ is real and
$\Gamma_{\rm im}=0,\,\Gamma_{\rm re}<0$, let $E_{1}=E_{4}$ for the Hermitian
model, for which $E_{2}>0,\,E_{3}>0$. We obtain that

$\Gamma_{\rm re}^{2}=4\Omega^{2}.$ (37)

The solution is $\Gamma_{\rm re}=-2\Omega=\Gamma_{c}$, which is the energy
level crossing point for the Hermitian model. Near critical point,
$E_{1}<E_{4}$ for $\Gamma_{\rm re}<\Gamma_{c}$, while $E_{4}<E_{1}$ for
$\Gamma_{\rm re}>\Gamma_{c}$. As fig.15 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$.

The order parameter $M^{x}$ of ground state can be given as

$M^{x}(\Gamma_{\rm
re}<\Gamma_{c})=\langle\psi_{R}^{1}|\sigma^{x}|\psi_{R}^{1}\rangle=0,$ (38)

and

$M^{x}(\Gamma_{\rm
re}>\Gamma_{c})=\langle\psi_{R}^{4}|\sigma^{x}|\psi_{R}^{4}\rangle=\frac{2\Omega(E+\Gamma_{\rm
re})}{2\Omega^{2}+(E+\Gamma_{\rm re})^{2}}.$ (39)

which tends to $2\Omega\Gamma_{c}/(2\Omega^{2}+\Gamma_{c}^{2})<0$ for
$\Gamma_{\rm re}\rightarrow\Gamma_{c}^{+}$. As shown in fig.16 of numerical
results, the number of spin-up and order parameter were rescaled with $L$ vary
with $\Gamma_{\rm re}$ for different system size. There are sudden jumps with
both observables near critical point. All results for different system size
collapse into same curve. We can conclude that order parameter $M^{x}$ of
ground state has sudden jump at two side of critical point, which implies
first-order phase transition.

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

$\hat{H}=\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}-\frac{1}{4}\Gamma_{\rm
re}\sum_{k=1}^{L}\hat{\sigma}_{z}^{k}+const.$ (40)

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$,
where agrees with critical point in Figs.15. 16 by numerics.