# Exploring the impact of fluctuation-induced criticality on non-hermitian
skin effect and quantum sensors

Clement Ehrhardt  Department of Physics, Stockholm University, AlbaNova
University Center, 106 91 Stockholm, Sweden  Jonas Larson  Department of
Physics, Stockholm University, AlbaNova University Center, 106 91 Stockholm,
Sweden

(August 28, 2024)

###### Abstract

In this paper, we present a concrete example that highlights how predictions
in non-Hermitian quantum mechanics can be inaccurately influenced by the
absence of environment-induced fluctuations in the model. Specifically, we
investigate the non-Hermitian skin effect and sensor in the Hatano-Nelson
model, contrasting it with a more precise Lindblad description. Our analysis
reveals that these phenomena can undergo breakdown when environmental
fluctuations come to the forefront, resulting in a non-equilibrium phase
transition from a localized skin phase to a delocalized phase. Beyond this
specific case study, we engage in a broader discussion regarding the
interpretations and implications of non-Hermitian quantum mechanics. This
examination serves to broaden our understanding of these phenomena and their
potential consequences.

###### pacs:

45.50.Pq, 03.65.Vf, 31.50Gh

##  I Introduction

In recent years, the field of non-Hermitian (NH) quantum mechanics (QM) has
experienced a remarkable resurgence ueda . This renaissance can be traced back
to the intriguing discovery that $\mathcal{PT}$-symmetric Hamiltonians, not
necessarily Hermitian, can yield real spectra bender . A pivotal moment in
this revival occurred with the introduction of biorthogonal QM bioqm , which
ignited debates about the fundamental nature of QM. It challenged the long-
held notion that observables must be represented solely by Hermitian operators
largeqm .

The focus of NH QM has evolved to explore novel phenomena that emerge when we
relax the constraints of Hermiticity and unitarity. One of the most
extensively studied phenomena is the NH skin effect skinref ; skineffect ;
skineffect2 , which renders extreme sensitivity to non-local perturbations ee
; expsens ; sensor . For certain NH local Hamiltonians with open boundary
conditions, all left/right eigenvectors $|\phi_{n}^{L,R}\rangle$ localize to
one of the edges, offering intriguing possibilities for detection of weak
signals sensor ; sensor4 .

NH QM often serves as an effective description of open quantum systems,
typically arising from the interaction with an external environment. However,
this approach raises questions about the treatment of fluctuations and the
potential violation of well-established quantum theorems nosignaling ;
nocloning ; LR ; entincrease ; unred .

This paper adopts a different perspective by employing the Lindblad master
equation (LME) as a foundational framework to analyze quantum systems exposed
to losses. Unlike NH QM, we do not neglect fluctuations, thus avoiding
concerns related to quantum jumps. We also explore the implications and
interpretations of NH theories in greater detail. Our study focuses on a
specific example, where fluctuations qualitatively alter the physics of the
system. We investigate a LME that reduces to the Hatano-Nelson (HN) model hn
in the absence of quantum jumps, revealing a breakdown of the NH skin effect
in favor of a delocalized phase. We discuss how such non-equilibrium
criticality relates to earlier models in the context of optical bistability
optbis ; optbis2 .

In conclusion, we examine the effects of fluctuations on NH QM, offering a
perspective that complements existing research comp1 ; physstate ; nhlindblad1
; nhlindblad2 ; nhlindblad2b ; hatanonelsonlindblad ; nhlindblad4 ; sensor2 ;
jumpperturb , especially by identifying a phenomenon of fluctuation-induced
criticality which qualitatively alters the physical properties. We aim to
provide a more detailed understanding of the role of fluctuations in NH QM.
This will help shed light on the applicabilities of the theory in the quantum
regime.

The paper is structured as follows: In the next section, we provide an in-
depth discussion of non-unitary time evolution, with a particular focus on its
description within the LME. We emphasize the importance of CPTP (Completely
Positive, Trace-Preserving) maps and use them to argue why eigenvectors of a
Liouvillian should not be considered as physical states. In Sec. III, we
introduce the model system, the HN model in Subsec. III.1, and its LME
realization in Subsec. III.2. Our main findings are presented in Sec. IV,
beginning with an exploration of the NH skin effect in Subsec. IV.1 and then a
discussion of how this translates to applications in sensing in Subsec. IV.2.
We conclude with a discussion in Sec. V. Additionally, we include two
appendices. The first provides general comments on open quantum systems
(appendix. A), and in the second, we demonstrate how our results for the HN
model also apply to the NH SSH model (appendix. B).

##  II Non-unitary quantum evolution

In this section, our primary aim is not to present new findings but rather to
provide context by discussing general aspects of quantum state evolution. We
specifically focus on the evolution generated by the LME. It is important to
note that the Liouvillian, which is responsible for governing time evolution
within the LME framework, is not represented by an observable. This
distinction leads to significant differences compared to Hamiltonian systems.
For instance, the eigenvectors of the Liouvillian do not typically represent
physical states.

Nonetheless, a NH Hamiltonian somewhat falls between a traditional Hamiltonian
and a master equation, but the field’s terminology tends to lean more towards
that of a Hamiltonian system. Having addressed these formal issues, the
subsequent section will explore a concrete example as we apply our knowledge
to the NH mode.

###  II.1 The Lindblad master equation

In Appendix A, we provide a more detailed description of open quantum systems.
In this section, we will simply state that our system, denoted as
$\mathcal{S}$, is weakly coupled to its surrounding environment. This
inevitably implies that the evolution of the system alone cannot be described
solely through unitary time evolution. However, we can still assume that the
state $\hat{\rho}(t)$, which characterizes the physical properties of the
system, adheres to the following physical state conditions:

$\begin{array}[]{lllll}(i)&&\mathrm{Tr}\left[\hat{\rho}(t)\right]=1,&&\mathrm{Normalization},\\\
\\\ (ii)&&\hat{\rho}(t)=\hat{\rho}^{\dagger}(t),&&\mathrm{Hermiticity},\\\ \\\
(iii)&&|\\!|\hat{\rho}(t)|\\!|\geq 0,&&\mathrm{Positivity}.\end{array}$ (1)

The first condition corresponds to standard normalization, which preserves
probabilities. The second condition ensures that all eigenvalues of
$\hat{\rho}(t)$ are real, while the third condition guarantees that all
eigenvalues are non-negative. These conditions are crucial for maintaining the
probability interpretation of quantum mechanics, e.g. avoidning negative
probabilities. It is important to note that these conditions must hold for all
times, $t$. A mapping of physical states into new physical states, hence
obeying the above conditions, is referred to as a completely positive trace-
preserving map (CPTP) nc ; preskill . In classical systems that emulate
quantum dynamics (for example after applying the paraxial approximation to
light propagation in non-linear media), deviations from the first condition
are analogous to the loss or gain of particles or intensity, leading to a
departure from the probabilistic nature.

Lindblad posed the question ingemar : What is the most general differential
equation describing the time evolution of a quantum state, ensuring that the
evolved state $\hat{\rho}(t)$ remains a valid density operator at all times?
The most general form of such a dynamical CPTP map (on a differential
structure) can be expressed in the Lindblad form lindblad

$\begin{array}[]{lll}\displaystyle{\frac{\partial}{\partial
t}}\hat{\rho}&=&\hat{\mathcal{L}}\left[\hat{\rho}\right]=\displaystyle{i\left[\hat{\rho},\hat{H}\right]+\hat{\mathcal{D}}\left[\hat{\rho}\right]}\\\
\\\
&=&\displaystyle{i\left[\hat{\rho},\hat{H}\right]+\sum_{k}\gamma_{k}\left(2\hat{L}_{k}\hat{\rho}\hat{L}_{k}^{\dagger}-\hat{L}_{k}^{\dagger}\hat{L}_{k}\hat{\rho}-\hat{\rho}\hat{L}_{k}^{\dagger}\hat{L}_{k}\right)}.\end{array}$
(2)

Here, we introduce the Liouvillian operator $\hat{\mathcal{L}}$, and the
dissipator operator $\hat{\mathcal{D}}$ accounts for the influence of the
environment. The $\gamma_{k}$ values represent the “decay rates” for channel
$k$, and the $\hat{L}_{k}$’s are the Lindblad jump operators bp . The above
LME can be put on the form

$\displaystyle{\frac{\partial}{\partial
t}}\hat{\rho}=\hat{\mathcal{L}}_{c}\left[\hat{\rho}\right]\equiv
i\left(\hat{\rho}\hat{H}_{\mathrm{eff}}-\hat{H}_{\mathrm{eff}}^{\dagger}\hat{\rho}\right)+\hat{\mathcal{J}}_{c}\left[\hat{\rho}\right],$
(3)

with the effective NH “Hamiltonian” defined as

$\hat{H}_{\mathrm{eff}}=\hat{H}-i\sum_{k}\gamma_{n}\hat{L}_{k}^{\dagger}\hat{L}_{k},$
(4)

and the jump super-operator

$\hat{\mathcal{J}}_{c}\left[\hat{\rho}\right]=2c\sum_{n}\gamma_{n}\hat{L}_{n}\hat{\rho}\hat{L}_{n}^{\dagger}.$
(5)

Please note that we use “Hamiltonian” in quotation marks because, in general,
it does not correspond to a traditional physical Hamiltonian. Specifically,
when using “Hamiltonian” we refer to Eq. (4). The subscript $0\leq c\leq 1$
parametrizes the master equation, where $c=1$ reproduces the correct LME (2),
and for $c=0$, the evolution is governed by the NH “Hamiltonian”
$\hat{H}_{\mathrm{eff}}$. We will exclusively focus on time-independent jump
operators. It is worth mentioning that a microscopic derivation of the LME (2)
typically relies on three approximations bp ; carmichael : the Markovian,
Born, and secular approximations. The properties of the environment and the
system-environment Hamiltonian determine the dissipator
$\hat{\mathcal{D}}\left[\hat{\rho}\right]$, including the rates $\gamma_{k}$
and the jump operators $\hat{L}_{k}$. Importantly, in especially cold atom and
optical systems, it is feasible, to a high degree, to engineer both the system
and its coupling to the environment to achieve desired Liouvillians lengineer
. In the following section, we provide explicit suggestions on how to utilize
relaxation to implement a HN-like model characterized by unbalanced left/right
hopping in a 1D tight-binding lattice.

###  II.2 Some general properties of the Lindblad master equation

An eigenvector, denoted as $\hat{\rho}_{j}$ and its corresponding eigenvalue,
$\mu_{j}$, of the LME are defined by the equation unik ; albert ; patrik

$\hat{\mathcal{L}}\left[\hat{\rho}_{j}\right]=\mu_{j}\hat{\rho}_{j}.$ (6)

In principle, we should refer to $\hat{\rho}_{j}$ as an “eigenmatrix”, but for
simplicity, we will continue to use the term “eigenvector”. This choice is
justified since the LME (2) can be vectorized as

$\frac{d}{dt}|\rho\rangle\\!\rangle=\hat{\mathcal{L}}_{v}|\rho\rangle\\!\rangle$
(7)

with the voctorized Liouvillian, $\hat{\mathcal{L}}_{v}$, defined as

$\begin{array}[]{lll}\displaystyle{\hat{\mathcal{L}}_{v}}&=&-i\left(\hat{H}\otimes\mathbb{I}+\mathbb{I}\otimes\hat{H}^{T}\right)\\\
\\\
&&+\displaystyle{\sum_{k}\gamma_{k}\\!\left[2\hat{L}_{k}\\!\otimes\hat{L}_{k}^{\dagger
T}-\hat{L}_{k}^{\dagger}\hat{L}_{k}\\!\otimes\mathbb{I}-\mathbb{I}\otimes\\!\left(\hat{L}_{k}^{\dagger}\hat{L}_{k}\right)^{T}\right]\\!.}\end{array}$
(8)

Now, given a finite Hilbert space dimension $D(\mathcal{H})=N$, $\hat{\rho}$
is represented as an $N^{2}$-component vector, $|\rho\rangle\\!\rangle$,
rather than an $N\times N$ matrix vec . The vector $|\rho\rangle\\!\rangle$
resides in a Liouville space $\mathcal{L}$, which is the direct product of two
Hilbert spaces, such that the dimension of the Liouville space is
$D(\mathcal{L})=N^{2}$.

For the numerical computations presented in Sec. IV, we utilize the vectorized
version of the LME. We select an appropriate (Fock) basis and express the
operators $\hat{H}$ and $\hat{L}_{k}$ as matrices, which are employed to
construct the Liouvillian matrix $\hat{\mathcal{L}}_{v}$. Subsequently, we
numerically diagonalize it to determine the Liouvillian spectrum $\mu_{j}$ and
corresponding eigenvectors
$|\rho_{j}\rangle\\!\rangle\leftrightarrow\hat{\rho}_{j}$ (including the
steady state). The steady state of Eq. (6) complies with
$\hat{\mathcal{L}}\left[\hat{\rho}_{\mathrm{ss}}\right]=0$. It is important to
note that the vector components of $|\rho\rangle\\!\rangle$ do not represent
probability amplitudes. For the scalar product, we have
$\langle\\!\langle\rho|\varrho\rangle\\!\rangle=\mathrm{Tr}\left[\hat{\rho}\hat{\varrho}\right]$.

The uniqueness of the steady state $\hat{\rho}_{\mathrm{ss}}$ has been
extensively discussed in previous works, starting with studies by Spohn and
continued by others; refer, for example, to Refs. unik . For LMEs possessing
symmetries, multiple steady states can arise, leading to non-trivial
situations reminiscent of spontaneous symmetry breaking, akin to continuous
phase transitions in closed systems albert ; patrik . We will return to this
in the next section when analyzing a LME-extension of the HN model.

The LME is a CPTP map, i.e. $dN/dt=0$ where
$N=\mathrm{Tr}\left[\hat{\rho}(t)\right]$, and
$\langle\psi|\hat{\rho}(t)|\psi\rangle\geq 0$ for any state $|\psi\rangle$.
However, it is important to note that as soon as the previously introduced
parameter $c\neq 1$, the CPTP property is generally lost. Specifically, in the
NH limit with $c=0$, it is only under very special circumstances that the
evolution remains trace-preserving, even when the spectrum is purely real. A
significant consequence of the CPTP property of the LME can be directly
inferred physstate ; sofia . First, we observe that an eigenvector evolves as

$\hat{\rho}_{j}(t)=e^{\mu_{j}t}\hat{\rho}_{j},$ (9)

and thus, its trace

$N_{j}(t)=\mathrm{Tr}\left[\hat{\rho}_{j}(t)\right]=e^{\mu_{j}t}\mathrm{Tr}\left[\hat{\rho}_{j}\right].$
(10)

However, since $dN(t)/dt=0$, we must have
$\mathrm{Tr}\left[\hat{\rho}_{j}\right]=0$ whenever $\mu_{j}\neq 0$. In other
words, every eigenvector, apart from the steady states, is traceless.
According to (1), this leads to the crucial result:

  1. Every eigenvector $\hat{\rho}_{j}$ of some Lindblad Liouvillian $\hat{\mathcal{L}}$, apart for its steady state(s), are unphysical.

For classical systems, it is well-established that only the steady state of a
master equation can exclusively consist of non-negative entries. In contrast,
all other eigenvectors possess at least one entry that is negative, rendering
them unsuitable for representing physical states, as these entries correspond
to probabilities vk . The result stated above parallels this observation, but
now applies to quantum systems. Consequently, this motivates us to express the
following:

  1. As with a master equation, the eigenvectors and eigenvalues of $\hat{\mathcal{L}}_{c}$ should not, in the general case, be interpreted as representing physical states and observable quantities, respectively.

This assertion is not confined solely to the LME case ($c=1$), as we contend
(as elaborated below) that the evolution generated by $\hat{H}_{\mathrm{eff}}$
should not be regarded as Hamiltonian time evolution. When one unravels the
Lindblad master equation in quantum trajectories unraveling , the NH
(renormalized) time evolution arises through post-selection. Thus, in this
case, while it may seem as if $\hat{H}_{\mathrm{eff}}$ generates the time
evolution, one should not forget that it is actually only true due to the
post-selection constraint.

In both classical and quantum scenarios, the eigenvectors still constitute a
complete set, apart at exceptional points, meaning that any other state can be
expressed as a linear combination of them,
$\hat{\rho}=\sum_{j}p_{j}\hat{\rho}_{j}$, using coefficients $p_{j}$. Thus,
the time-evolved state can be expressed as

$\hat{\rho}(t)=\hat{\rho}_{\mathrm{ss}}+\sum_{j}p_{j}e^{\mu_{j}t}\hat{\rho}_{j},$
(11)

where it is implied that the sum excludes eigenvectors with $\mu_{j}=0$. To
ensure CPTP behavior, the eigenvalues must satisfy $\mathrm{Re}(\mu_{j})\leq
0$ bp ; huelga ; nonneg . Consequently, all states $\hat{\rho}_{j}$ with a
non-zero real part $\mathrm{Re}(\mu_{j})$ are exponentially suppressed as time
progresses. The Liouvillian gap, defined as lgap

$\Delta_{\mathcal{L}}=\min_{j}\mathrm{Re}(\mu_{j}),$ (12)

can be regarded as providing an initial estimate for the relaxation time scale
towards the steady state (although the scenario can be more intricate lgap2 ).

In the context of quantum phase transitions, in the thermodynamic limit, the
ground state exhibits non-analytic behavior at the critical point, and
additionally, the spectrum necessarily becomes gapless at this juncture.
Analogously, the steady state $\hat{\rho}_{\mathrm{ss}}$ may display similar
non-analytic behavior, accompanied by the vanishing of the Liouvillian gap
lgap ; lcrit . It is important to note that the steady state does not decay,
although under certain conditions nons , $\mathrm{Re}(\mu_{j})$ may equal 0
while $\mathrm{Im}(\mu_{j})\neq 0$, potentially resulting in non-stationary
states. However, this requires that the eigenvalues $\mu_{j}$ must appear in
complex conjugate pairs if they possess a non-vanishing imaginary part.

Returning to Eq. (3), if we put the last term to zero (i.e. $c=0$) the
evolution is governed by the NH “Hamiltonian” $\hat{H}_{\mathrm{eff}}$ of Eq.
(4). In this case, let us assume that we know the right eigenvectors

$\hat{H}_{\mathrm{eff}}^{\dagger}|\varphi_{l}^{R}\rangle=\nu_{n}|\varphi_{l}^{R}\rangle.$
(13)

That is, the left eigenvectors obey

$\hat{H}_{\mathrm{eff}}|\varphi_{l}^{L}\rangle=\nu_{n}^{*}|\varphi_{l}^{L}\rangle.$
(14)

The eigenvalues and eigenvectors of $\hat{\mathcal{L}}_{c=0}$ become

$\mu_{j}^{0}=i(\nu_{l}^{*}-\nu_{k}),\hskip
28.45274pt\hat{\rho}_{j}^{0}=|\varphi_{k}^{R}\rangle\langle\varphi_{l}^{R}|.$
(15)

Here, the superscript 0 denotes the case $c=0$, and $j$ replaces the double
indices $(l,k)$. As mentioned earlier, if the Hilbert space dimension is
finite, $D(\mathcal{H})=N$, the Liouville space dimension is $N^{2}$, as seen
in Eq. (15) since $1\leq k,\,l\leq N$ implies $1\leq j\leq N^{2}$. Thus, we
have far more eigenvectors/values of the Liouvillian $\hat{\mathcal{L}}_{c=0}$
than for the NH “Hamiltonian” $\hat{H}_{\mathrm{eff}}$.

We can parametrize the eigenvectors of $\hat{\mathcal{L}}_{c}$ with the number
$c$ (i.e. $\hat{\rho}_{j}^{c}$), such that $\hat{\rho}_{j}^{c=0}$ reproduces
the eigenvectors in (15), while $\hat{\rho}_{j}^{c=1}$ provides those of Eq.
(6). For finite dimensionsional Hilbert spaces we expect the vectors
$\hat{\rho}_{j}^{c}$ and eigenvalues $\mu_{j}^{c}$ to be analytic in $c$, with
the exceptions at possible exceptional points ep . Letting $c=0$, if $\nu_{l}$
is real, it follows that the vectors
$\hat{\rho}_{j}^{0}=|\varphi_{l}^{R}\rangle\langle\varphi_{l}^{R}|$ are steady
states and also physical. However, if $\mathrm{Im}(\nu_{j})\neq 0$, these
eigenvectors
$\hat{\rho}_{j}^{0}=|\varphi_{l}^{R}\rangle\langle\varphi_{l}^{R}|$ are no
longer steady states since their norms are not preserved.

Let us provide some context for the discussion above. As mentioned in the
introduction, NH QM challenges well-established physical concepts.
Furthermore, an ongoing debate surrounds the clear interpretation of the
theory emerging from NH QM, particularly with regard to which states should be
employed to describe time evolution hatanonelsonlindblad ; rrdebate1 .

One approach to circumvent potential issues, such as violating the no-
signaling theorem, is the introduction of biorthogonal QM nogo ; bioqm . This
approach hinges on the biorthogonality property, which allows for the
construction of mutually orthogonal left $|\varphi_{j}^{L}\rangle$ and right
$|\varphi_{j}^{R}\rangle$ eigenvectors, such that
$\langle\varphi_{l}^{L}|\varphi_{j}^{R}\rangle=\delta_{lj}$. For example, this
property leads to the modified resolution of identity, which becomes
$\hat{O}=\sum_{j}o_{j}|\phi_{j}^{R}\rangle\langle\phi_{j}^{L}|$. Additionally,
the spectral resolution of an operator can be expressed as
$\hat{O}=\sum_{j}o_{j}|\phi_{j}^{R}\rangle\langle\phi_{j}^{L}|$, where $o_{j}$
and $|\phi_{j}^{R}\rangle$ ($\langle\phi_{j}^{L}|$) represent the eigenvalue
and right (left) eigenvectors of the operator, respectively. In the
biorthogonal framework, assuming an initial state $|\psi(0)\rangle$, the state
at a later time is described by two vectors: $|\psi^{L,R}(t)\rangle$, referred
to as the ‘left’ and the ‘right’ evolved state. Expectations of observables
$\hat{O}$ should also be evaluated according to this ‘state,’ given by

$\mathcal{O}_{|\psi^{L,R}\rangle}(t)=\langle\psi^{L}(t)|\hat{O}|\psi^{R}(t)\rangle=\mathrm{Tr}\left[\hat{O}\hat{\rho}(t)\right],$
(16)

with $\hat{\rho}(t)=|\psi^{R}(t)\rangle\langle\psi^{L}(t)|$.

This can be extended further; the expectation of an operator for the $j$’th
eigenvector becomes
$\mathcal{O}_{j}=\langle\varphi_{j}^{L}|\hat{O}|\varphi_{j}^{R}\rangle$. When
applied to the position operator $\hat{n}_{n}=|n\rangle\langle n|$ (where
$|n\rangle$ represents the particle localized to site $n$) of the HN model,
one finds that the “biorthogonal” eigenvectors are not localized at the edges
but instead are delocalized within the bulk el2  (see Eq. (24) below for the
eigenvectors).

Of course, if the generator of time evolution is a Hermitian Hamiltonian, this
reproduces standard quantum mechanics since $\langle\psi^{L}|$ will equal the
right bra vector $\langle\psi^{R}|$. However, while the biorthogonal approach
seems to resolve some issues, it comes with caveats. First, we note that
$\langle\varphi_{l}^{L}|\varphi_{j}^{R}\rangle=\delta_{lj}$ does not set the
norm of the left/right eigenvectors independently of each other. In fact, this
results in a family of scalar products parametrized by a metric $\eta$ nogo ;
elisabet . Thus, an ambiguity, similar to a gauge freedom, arises. It has been
argued that for a given Hamiltonian, there is a preferable metric to be used,
which, however, implies that the metric, and thereby the scalar product,
changes when, for example, you add one particle or site (modifying the Hilbert
space dimension) to your system bioqm ; elisabet .

Secondly, the fact that a physical state must be ascribed both a bra- and ket-
vector is non-intuitive using common knowledge, and we have not even addressed
mixed states. Now, it has been argued that the biorthogonal scalar product is
not the one to be used in order to describe the time evolution of actual
physical systems rrdebate1 ; hatanonelsonlindblad . This aligns with our
description. If we take the LME as our starting point and think that we can,
at least in theory, connect it to NH QM by ‘turning off’ the jump terms, we
should not end up with the biorthogonal left/right formalism, but rather with
a right/right (or left/left) density operator com1 . In the post-selection
approach, it is yet again the right/right state (complemented with
renormalization) that describes the system.

It turns out that the spectral properties (15) of the NH Hamiltonian
$\hat{H}_{\mathrm{eff}}$ can capture the full Liouvillian spectrum. To be more
precise, under specific conditions, the spectrum of $\hat{\mathcal{L}}_{c=0}$
can be directly mapped to the spectrum of $\hat{\mathcal{L}}_{c=1}$ lspec .
This mapping becomes feasible when the system Hamiltonian supports particle
conservation, as indicated by $\left[\hat{N},\hat{H}\right]=0$, and every jump
operator adheres $\left[\hat{L}_{k},\hat{N}\right]=\hat{L}_{k}$. This
situation is relevant in the context of spontaneous decay bp . It is worth
noting that our Lindblad representation of the HN model, as presented in the
following section, does not satisfy this second condition, and therefore, we
cannot utilize such a property.

Another interesting scenario arises when the Liouvillian has a quadratic
dependence on the creation/annihilation operators
$\hat{a}_{k}^{\dagger}/\hat{a}_{k}$ prosen . Similar to quadratic
Hamiltonians, the quadratic form of the Liouvillian can be employed for
diagonalization through a generalized Bogoliubov–de Gennes approach. This
method has been recently applied in the investigation of the skin effect in
Liouvillians nhlindblad2b . However, it is important to note that the physical
realization of the HN model we have in mind, as described in Subsec. III.2,
does not fall within this class of solvable models.

##  III Model system: Liouvillian for the Hatano-Nelson model

###  III.1 The Hatano-Nelson model

A fundamental model often used in the study of NH QM is the one proposed by
Hatano and Nelson. This model describes a single particle within a one-
dimensional tight-binding lattice of N sites, where the hopping between sites
is unbalanced in the left and right directions hn ; hn2 . For open boundary
conditions (BC), the model is represented by the (NH) Hamiltonian as

$\hat{H}_{\mathrm{HN}}=\sum_{n=1}^{N-1}\left[\left(1-\delta\right)\hat{a}_{n}^{\dagger}\hat{a}_{n+1}+(1+\delta)\hat{a}_{n+1}^{\dagger}\hat{a}_{n}\right].$
(17)

Here, $\hat{a}_{n}$ ($\hat{a}_{n}^{\dagger}$) represent the annihilation
(creation) operators for a particle at site $n$, satisfying the single
particle constraint
$\hat{N}=\sum_{n=1}^{N}\hat{a}_{n}^{\dagger}\hat{a}_{n}\equiv\sum_{n=1}^{N}\hat{n}_{n}=1$.
The parameter $\delta$ varies between 0 and 1 and indicates the degree of
asymmetry in the hopping to the left and right.

For a single particle, the model can also be expressed using bracket notation,
for example,
$\hat{a}_{n+1}^{\dagger}\hat{a}_{n}\leftrightarrow|n+1\rangle\langle n|$, and
so forth, where $|n\rangle$ represents the number state of the particle
localized at the $n$-th site. The model is typically presented in
dimensionless units, scaled with respect to the symmetric hopping amplitude.
It can be convenient to decompose the Hamiltonian into ‘real’ and ‘imaginary’
parts as follows

$\hat{H}_{\mathrm{HN}}=\hat{H}_{R}+i\delta\hat{H}_{I},$ (18)

where

$\hat{H}_{R}=\sum_{n=1}^{N-1}\left(\hat{a}_{n}^{\dagger}\hat{a}_{n+1}+\hat{a}_{n+1}^{\dagger}\hat{a}_{n}\right)$
(19)

and

$\hat{H}_{I}=i\sum_{n=1}^{N-1}\left(\hat{a}_{n}^{\dagger}\hat{a}_{n+1}-\hat{a}_{n+1}^{\dagger}\hat{a}_{n}\right)$
(20)

Both $\hat{H}_{R}$ and $\hat{H}_{I}$ are hermitian. Notably, in the case of
open boundary conditions, the commutator
$\left[\hat{H}_{R},\hat{H}_{I}\right]$ vanishes, except for the first and last
diagonal elements.

For future reference, it is informative to introduce the ladder operators,
denoted as $\hat{E}=\sum_{n}\hat{a}_{n}^{\dagger}\hat{a}_{n+1}$ and its
hermitian conjugate $\hat{E}^{\dagger}$. These operators satisfy the following
relationships: $\hat{E}|n\rangle=|n-1\rangle$ and
$\hat{E}^{\dagger}|n\rangle=|n+1\rangle$ (it should be noted that for a finite
lattice, $\hat{E}|1\rangle=0$ and $\hat{E}^{\dagger}|N\rangle=0$). When
combined with the operator $\hat{E}_{0}$, which follows the relation
$\hat{E}_{0}|n\rangle=n|n\rangle$, these three operators collectively
constitute what is known as the Euclidean algebra euclid

$\left[\hat{E},\hat{E}_{0}\right]=-\hat{E},\hskip
14.22636pt\left[\hat{E}^{\dagger},\hat{E}_{0}\right]=\hat{E}^{\dagger},\hskip
14.22636pt\left[\hat{E},\hat{E}^{\dagger}\right]=0.$ (21)

Alternatively, we can introduce the “psoition” and “momentum” operators
$\hat{E}_{x}=\hat{E}+\hat{E}^{\dagger}$ and
$\hat{E}_{p}=i\left(\hat{E}-\hat{E}^{\dagger}\right)$, and the “Hamiltonian”
reads $\hat{H}_{\mathrm{HN}}=\hat{E}_{x}+i\delta\hat{E}_{p}$.

For periodic BC, the real and imaginary parts commute and the eigenvectors are
delocalized bulk states. For open BC the spectrum reads ee ; hnspec

$\nu_{j}=2\sqrt{1-\delta^{2}}\cos(j\pi/(N+1)),$ (22)

while for periodic BC one has

$\nu_{j}=2\left[\cos(2\pi j/N)-i\delta\sin(2\pi j/N)\right],$ (23)

where $j=1,2,\dots,N$. Consequently, the spectrum is purely real for open
boundary conditions. However, the right and left eigenvectors are not
orthogonal, which leads to non-unitary time evolution. All right eigenvectors,
while not normalized, are exponentially localized towards the right edge and
can be represented in terms of the number states as:

$|\varphi_{j}^{R}\rangle=\sum_{n=1}^{N}\left(\frac{1+\delta}{1-\delta}\right)^{n}\sin\left(\frac{nj\pi}{N+1}\right)|n\rangle.$
(24)

The left eigenvectors can be obtained by substituting $\delta$ with $-\delta$.
By using this substitution, we find that the biorthogonal QM provides the
eigenvector site occupations
$P_{j}(n)=\langle\varphi_{j}^{L}|\hat{n}_{n}|\varphi_{j}^{R}\rangle\propto\sin^{2}\left(\frac{nj\pi}{N+1}\right)$.
The bulk nature of these states, when combining left and right eigenvectors,
was previously observed by Hatano and Nelson hn2 .

###  III.2 Lindblad implementation of the Hatano-Nelson model

Our objective is to construct jump operators $\hat{L}_{k}$, which enable the
NH “Hamiltonian” in Eq. (4) to match the HN Hamiltonian (17). Once these
operators are identified, we propose a physical system where this can be
implemented. By separating the real and imaginary components of both
equations, we should obtain the following relationships

$\hat{H}_{R}=\hat{H}$ (25)

and

$\delta\hat{H}_{I}=-\sum_{k}\gamma_{k}\hat{L}_{k}^{\dagger}\hat{L}_{k}.$ (26)

The first equation implies that our Hamiltonian should resemble a simple
$N$-site tight-binding chain. The second identity is more complex, as the
spectrum of $\hat{H}_{I}$ is not strictly non-negative, while
$\hat{L}_{k}^{\dagger}\hat{L}_{k}$ is positive semi-definite (recall
$\gamma_{k}\geq 0$). Nonetheless, we can resolve this issue by shifting the
entire spectrum of the NH “Hamiltonian” by an imaginary constant of $2i\delta$
hatanonelsonlindblad . It is worth noting that the decomposition into jump
operators is not unique; instead, there are numerous possible combinations. To
illustrate this, the right-hand side of equation (26) can be compactly
expressed as $-\hat{\bf L}^{\dagger}\Gamma\hat{\bf L}$, where $\hat{\bf
L}=(\hat{L}_{1},\,L_{2},\dots)^{t}$ is a column vector containing all jump
operators (for a finite dimension $D(\mathcal{H})=N$, the number of
independent jump operators can be limited to $N^{2}-1$), and $\Gamma$ is a
diagonal matrix with $\gamma_{k}$ along its diagonal. It is evident that this
term remains invariant under some unitary transformation $\hat{U}$, i.e.,
$\hat{\bf L}^{\prime}=\hat{U}\hat{\bf L}$ and
$\Gamma^{\prime}=\hat{U}\Gamma\hat{U}^{-1}$. The specific choice should be
determined by the physical system in question. Two cases, in particular, come
to mind:

  1. 1.

Local decay channels. At each site $n$, we can associate a jump operator
$\hat{L}_{n}$. When this operator acts on the state $|n\rangle$, there is a
non-negligible probability for the particle to transition to $|n+1\rangle$. A
straightforward approach might be to set
$\hat{L}_{n}=\hat{a}_{n}\hat{a}_{n+1}^{\dagger}$. However, this does not
suffice, as $\hat{L}_{n}^{\dagger}\hat{L}_{n}$ results in $\hat{n}_{n}$, which
cannot be used to construct the desired effective NH “Hamiltonian”. Instead,
we should employ the following expression:

$\hat{L}_{n}=i\hat{a}_{n}\hat{a}_{n+1}^{\dagger}+\hat{n}_{n},$ (27)

ensuring that the application of the jump operator to $|n\rangle$ leaves the
state in a superposition of $|n\rangle$ and $|n+1\rangle$.

  2. 2.

Collective decay channels. If the jump operator does not “keep track” of the
particle’s position, we say that the jump occurs independently of the
particle’s position. In this case, we sum the local jump operators into a
single one, i.e.

$\hat{L}_{c}=i\sum_{n=1}^{N-1}\hat{a}_{n}\hat{a}_{n+1}^{\dagger}+\sum_{n=1}^{N}\hat{n}_{n}\equiv\hat{E}+\mathbb{I},$
(28)

where $\hat{E}$ was defined above Eq. (21).

The two scenarios described above result in different physical realizations.
Unraveling the LME in terms of stochastic ’quantum trajectories’ unraveling
provides a physical picture of the system’s evolution. In summary, when we
possess complete knowledge about the environment, the system should evolve
deterministically with the NH “Hamiltonian” (4). However, this evolution is
occasionally interrupted by quantum jumps, whose effects are described by the
application of jump operators to the state. This process involves the
instantaneous renormalization of the state under non-unitary time evolution.

In practice, ‘keeping track’ of the environment involves monitoring any
photons spontaneously emitted by the system. If a photon is detected, we can
infer that the system has undergone a stochastic jump. This provides a
distinct physical differentiation between the two scenarios. In the latter
situation, the emitted photon does not yield information about the particle’s
position. In contrast, in the former scenario, a recorded photon not only
signifies that a jump has occurred but also reveals the particle’s location,
effectively collapsing the wave function within the lattice.

As a result, it becomes evident that we can anticipate different behaviors
under the evolution driven by the different Lindblad operators. For the
remainder of our discussion, we will focus on the collective decay channel.
However, for a finite lattice comprising $N$ sites, we should note that the
jump operator $\hat{L}_{c}$ of Eq. (28) is insufficient. We need to supplement
it with an additional dephasing operator
$\hat{L}_{1}=\hat{a}_{1}^{\dagger}\hat{a}_{1}$ for the first site. With this
consideration in mind we have

$\hat{L}_{c}^{\dagger}\hat{L}_{c}+\hat{L}_{1}^{\dagger}\hat{L}_{1}=\hat{H}_{I}+2\mathbb{I},$
(29)

and by further identifying $\gamma_{c}=\gamma_{1}\equiv\gamma=-\delta$ we have
our LME corresponding to the HN model

$\begin{array}[]{lll}\displaystyle{\frac{\partial}{\partial
t}}\hat{\rho}&=&\displaystyle{i\left[\hat{\rho},\hat{H}\right]+\gamma\left(2\hat{L}_{c}\hat{\rho}\hat{L}_{c}^{\dagger}-\hat{L}_{c}^{\dagger}\hat{L}_{c}\hat{\rho}-\hat{\rho}\hat{L}_{c}^{\dagger}\hat{L}_{c}\right)}\\\
\\\
&&+\displaystyle{\gamma\left(2\hat{L}_{1}\hat{\rho}\hat{L}_{1}^{\dagger}-\hat{L}_{1}^{\dagger}\hat{L}_{1}\hat{\rho}-\hat{\rho}\hat{L}_{1}^{\dagger}\hat{L}_{1}\right)},\end{array}$
(30)

where the Hamiltonian $\hat{H}$ identifies the real part (18) of the HN
Hamiltonian.

Figure 1: Illustration of the Implementation of the HN Model as an open
quantum system. In this schematic representation, atomic lattice sites are
depicted as black and red dots, signifying two distinct internal hyperfine
levels $|g\rangle$ and $|e\rangle$. These atoms are confined within an optical
lattice, which is state-dependent. The lattice is designed so that the two
internal states experience slight shifts concerning each other. Additionally,
the excited atomic state $|e\rangle$ is inherently unstable and quickly
relaxes (with an atomic rate $\gamma_{a}$) back to the stable atomic state
$|g\rangle$. As a result of this shift in the respective potential minima and
the relaxation process, an effective drift occurs in the lattice, moving in
the right direction. Specifically, when an atom in the excited state
$|e\rangle$ at site $n$ , it relaxes back to the ground state $|g\rangle$ at
either the site $n$ or $n+1$, while relaxation to other sites is negligible.

In this section, we offer a conceptual overview of how this model could be
implemented within a cold atom setup coldatom . A similar concept was recently
proposed using motional sidebands in a trapped ion setup trapsetup . Our
discussion centers on an atom that possesses two internal hyperfine levels:
$|g\rangle$ (ground) and $|e\rangle$ (excited). This atom is confined within a
one-dimensional lattice. Under the usual approximations, which include tight-
binding and single-band assumptions, we implement a resonant classical drive
between the two internal atomic states. This leads to the following second-
quantized lattice Hamiltonian

$\begin{array}[]{lll}\hat{H}&=&\displaystyle{-t_{g}\sum_{n}\left(\hat{a}_{n}^{\dagger}\hat{a}_{n+1}+h.c.\right)-t_{e}\sum_{n}\left(\hat{b}_{n}^{\dagger}\hat{b}_{n+1}+h.c.\right)}\\\
\\\
&&\displaystyle{+g\sum_{n}\left(\hat{b}_{n}^{\dagger}\hat{a}_{n}+h.c.\right)}.\end{array}$
(31)

Here, $t_{g}$ and $t_{e}$ represent the tunneling rates of the two atomic
species, and $\hat{a}_{n}$ ($\hat{a}_{n}^{\dagger}$) and $\hat{b}_{n}$
($\hat{b}_{n}^{\dagger}$) correspond to the single-site annihilation
(creation) operators. The parameter $g$ signifies the effective Rabi coupling.
Notably, the driving mechanism couples atomic internal states within a single
site. This implies a smooth laser profile and a sufficiently deep lattice. For
simplicity, we can assume $|t_{e}|\ll|t_{g}|$ to suppress lattice dynamics of
the excited states. We normalize energies in terms of the tunneling rate
$t_{g}$, setting $t_{g}=1$ from this point forward. If the $|e\rangle$-level
rapidly relaxes to the ground state $|g\rangle$, adiabatic elimination is
applicable adel . Consequently, we commence with the above Hamiltonian,
complemented by a bath of oscillators inducing couplings between the two
atomic levels. In this derivation, we employ the standard approximations,
which include Born, Markov, and secular approximations bp . The physical setup
is depicted in Figure 1. In the resulting LME for the atomic ground state
$|g\rangle$, the effective “decay rate” $\gamma$ is directly proportional to
the Rabi frequency $g$. It also significantly depends on the Franck-Condon
factors, which are proportional to the overlaps between Wannier functions of
the two species. The alignment of the two lattices and their lattice depths
should ensure that a Wannier function localized to the $n$’th site of the
“red” lattice primarily overlaps with Wannier functions at the $n$’th and
$n+1$’th lattice sites of the “black” lattice. In this scenario, and under the
assumption of precise phase alignment between the two decay channels (see Eq.
(28)), we achieve the desired LME.

##  IV Case studies

In the section we consider two phenomena within the HN model: the skin effect
skinref ; skineffect ; skineffect2  and NH sensors sensor ; sensor2 ; sensor3
. The aim is to compare and contrast outcomes obtained using either the HN
model of Eq. (17) or the full LME (30).

Figure 2: (Color online) In the upper plot (a) we display the scaled position
(32) for both the HN model (black), calculated as
$\xi=\langle\varphi_{N}^{R}|\hat{S}|\varphi_{N}^{R}\rangle/N$, and the full
LME (red). The solid lines represent a lattice of length $N=41$, while the
dashed lines represent $N=21$. In the lower plot (b), we illustrate the purity
(33) of the steady state $\hat{\rho}_{\mathrm{ss}}$ (purple) and the fidelity
(34) between the HN “ground state” and the Lindblad steady states (blue). The
key observation is the appearance of a critical point at $\gamma_{c}=1/2$ in
the LME, which is absent in the HN model. Below the critical point,
$\gamma<\gamma_{c}$, fluctuations become essential, leading to the
destabilization of the skin effect and system delocalization. Above the
critical point,$\gamma>\gamma_{c}$, the skin effect persists, even though the
fidelity is not very high.

###  IV.1 Skin effect

While the phenomenon was initially discovered by Hatano and Nelson in the
1990s hn ; hn2 , the term “non-hermitian skin effect” was coined by Yau and
Wang more than two decades later when they conducted a more in-depth and
comprehensive study skinref . When considering a NH lattice “Hamiltonian” with
open boundary conditions, the skin effect is characterized by the eigenvectors
becoming localized at the edges. This has a profound consequence: the system
becomes exponentially sensitive to non-local perturbations hn2 . This
sensitivity can be demonstrated by extending standard perturbation theory to
NH models sensor . Hence, NH edge localization and exponential sensitivity are
two sides of the same coin. We will revisit the latter concept in Subsec.
IV.2.

To characterize edge localization, we consider a lattice with $N=2n-1$ (where
$n$ is a positive integer) sites and introduce the matrix
$\hat{S}=\mathrm{diag}(-n:n)$ along with the scaled position expectation $\xi$
defined as

$\xi=\frac{1}{N}\mathrm{Tr}\left[\hat{\rho}_{\mathrm{ss}}\hat{S}\right].$ (32)

If $\xi=\pm 1$, the state $\hat{\rho}$ is maximally localized at one of the
edges, whereas $\xi=0$ indicates that the state is centered within the
lattice. The uncertainty
$\Delta\xi=\sqrt{\langle\hat{S}^{2}\rangle-\langle\hat{S}\rangle^{2}}$
determines the degree of localization of the state. Therefore, when
$\Delta\xi$ is of the order of $N$, the state becomes delocalized over the
entire lattice. The purity of the state is quantified by

$P=\mathrm{Tr}\left[\hat{\rho}_{\mathrm{ss}}^{2}\right].$ (33)

Notably, in the context of NH QM, the (normalized) eigenvectors are pure, and
their normalization ensures that $P_{\mathrm{NH}}=1$.

In the thermodynamic limit of an infinite lattice, where
$\left[\hat{L}_{c},\hat{L}_{c}^{\dagger}\right]=0$, the steady state of the
LME simplifies to the maximally mixed state,
$\hat{\rho}_{\mathrm{ss}}\propto\mathbb{I}$ patrik . However, for a finite
lattice, non-trivial steady states can emerge, particularly those localized at
the edges. When comparing this with the HN model, the question arises
regarding which state serves as the counterpart of the steady state in the HN
model.

One potential approach is to consider the right eigenvector
$|\varphi_{j}^{R}\rangle$ with the largest imaginary part of its eigenvalue,
i.e., $\max\left[\mathrm{Im}(\nu_{j})\right]$. This vector would represent the
steady state of $\hat{H}_{\mathrm{HN}}$, provided that we renormalize it under
time-evolution. However, for open boundary conditions, the spectrum is real.
Instead, we consider the eigenvector with the smallest eigenvalue, which
corresponds to the state $|\varphi_{N}^{R}\rangle$ as defined in Eq. (24). In
this way, we argue that this state, in some sense, mimics the ground state of
the system. The fidelity between this state and the steady state of the LME is
defined as

$F=\langle\varphi_{N}^{R}|\hat{\rho}_{\mathrm{ss}}|\varphi_{N}^{R}\rangle,$
(34)

The spectrum and eigenvectors of the LME are determined through exact
diagonalization of the vectorized Liouvillian (7). The numerical results for
the three defined quantities are presented in Fig. 2. In the upper plot (a),
we display the position (32) for both the LME and the HN model. We vary the
rate $\gamma$ and consider two system sizes. In the lower plot (b), we provide
the purity (33) of the steady state and the fidelity (34). There is a
significant distinction between the two models. The HN model exhibits edge
localization, or the skin effect, for all $\gamma$ values, and for $\gamma=1$,
only the rightmost site is populated.

In contrast, the LME supports two different phases: a delocalized phase for
$\gamma<\gamma_{c}=1/2$ and a localized phase for $\gamma>\gamma_{c}$. In the
localized phase, the state primarily occupies the right edge, resembling a
skin state, while in the delocalized phase, the state spreads across the
entire lattice. The purity, as seen in (b), illustrates that the delocalized
phase is approximately a maximally mixed state. Although the states
$\hat{\rho}_{\mathrm{ss}}$ and $|\varphi_{N}^{R}\rangle$ exhibit similar
localization properties near $\gamma=1$, the fidelity reveals that they are,
in fact, quite distinct.

Even though only the steady state represents a physical state, we can still
analyze the properties of the remaining eigenvectors of the Liouvillian. Not
shown here, but for $\gamma>\gamma_{c}$, they tend to localize at the edge,
while for $\gamma<\gamma_{c}$, they become delocalized.

To provide insight into the critical properties of the delocalized-to-
localized phase transition, we plot the real part of the Liouvillian spectrum
$\mu_{j}$ in Fig. 3 (a). It is evident that the Liouvillian gap
$\Delta_{\mathcal{L}}$, as defined in Eq. (12), closes precisely at the
critical point $\gamma_{c}$. In this figure, we consider a lattice with $N=41$
sites, resulting in a Liouvillian matrix of dimensions $1681\times 1681$. In
the thermodynamic limit, the real part of the spectrum becomes gapless and
continuous within the delocalized phase.

For a critical model, the closure of the gap concerning the system size
typically follows universal scaling:

$\Delta_{\mathcal{L}}\sim N^{-1/\nu},$ (35)

where $\nu$ represents the correlation length critical exponent cardy . Our
numerical findings affirm that the exponent $\nu$ is 1/2, as demonstrated in
Fig. 2 (b). This value corresponds to the exponent of a mean-field critical
model, suggesting that quantum fluctuations are suppressed in comparison to
fluctuations arising from the environment.

It might appear counterintuitive that the “fluctuation-induced breakdown” of
the skin effect occurs for weak couplings ($\gamma$) rather than when the
system is strongly coupled to its environment. This can be understood by
considering the steady state, which describes the system after an infinitely
long time. In the delocalized phase, the system is gapless, and the relaxation
time diverges. Consequently, there is more room for fluctuations to manifest
compared to the localized phase. This phenomenon resonates with what is known
in adiabatic quantum computing and quantum control in atomic physics stirap .

Figure 3: Real part of the Liouvillian spectrum $\mu_{j}$ (a) as a function of
$\gamma$ for a fixed system size $N=41$ (a). The lower plot (b) demonstrate
the gap closening ($\Delta_{\mathcal{L}}\rightarrow 0$ as
$N\rightarrow\infty$) at the critical point $\gamma=1/2$. It is evident from
the upper plot how the Liouvillian gap, $\Delta_{\mathcal{L}}$, opens up
beyond $\gamma=1/2$, while for $\gamma<1/2$, the numerical results indicates a
gapless, continuous spectrum. From (b), we find the power-law scaling with an
exponent $\nu=1/2$ characteristic of a for a mean-field critical point.

For reasons that will become clear, let us introduce a generalized Fock space.
The Liouvillian Fock states, denoted as $|l\rangle\\!\rangle$ and residing in
the Liouvillian space, can be identified through vectorization as

$|n\rangle\langle m|\rightarrow|N(n-1)+m\rangle\\!\rangle,$ (36)

such that $1\leq l\leq N^{2}$. The vectorized LME is provided in (7), and its
formal solution is represented as

$|\rho(t)\rangle\\!\rangle=\exp\left(\hat{\mathcal{L}}_{v}t\right)|\rho(0)\rangle\\!\rangle$
(37)

. We can now extend the concept of Fock state lattices pil  to Liouvillian
Fock state lattices. The idea here is to envision the Liouvillian matrix
$\hat{\mathcal{L}}_{v}$ as describing hopping in a lattice, with its sites
representing the Fock states $|l\rangle\\!\rangle$. In other words, the
components of the vector $|\rho\rangle\\!\rangle$ correspond to the
populations at different lattice sites. Specifically, the diagonal elements
$\varepsilon_{i}$, defined as $\varepsilon_{i}\equiv\langle\\!\langle
l|\hat{\mathcal{L}}_{v}|l\rangle\\!\rangle$, represent on-site “energies”,
while the off-diagonal elements $\langle\\!\langle
l|\hat{\mathcal{L}}_{v}|k\rangle\\!\rangle$ characterize the tunneling
amplitudes between sites $l$ and $k$. Importantly, $\varepsilon_{l}$ is real
and $\varepsilon_{l}\leq 0$, implying that the diagonal terms induce on-site
dissipation.

In the context of the problem at hand, it turns out that the Liouvillian Fock
state lattice takes the form of a square lattice, as illustrated in Fig. 4.
The thick dots represent the lattice sites, and the arrows indicate the
allowed tunnelings between these sites. The color shading of the sites
reflects the magnitude of on-site dissipation, with a decreasing order from
black to gray to white.

Having identified the Fock state lattice, we make the following observations.
There is no asymmetry in the horizontal/vertical tunnelings, as one might
expect from an HN model. However, there is diagonal tunneling that is non-zero
only in one direction. In addition to the imaginary on-site terms, these
diagonal directional tunneling terms imply that the Liouvillian matrix becomes
NH. Furthermore, these diagonal processes result from quantum jumps, driven by
fluctuations, and consequently, they do not appear in the lattice emerging
from the HN model. More precisely, the corresponding HN lattice is simply a
two-dimensional version of the HN model with imbalanced left/right tunnelings.
Hence, quantum jumps not only induce diagonal tunneling terms but also alter
the nearest neighbor tunneling amplitudes, making them balanced.

As previously mentioned, the on-site “energies” vary in the lattice. For
instance, at site $l=N^{2}$ (corresponding to the rightmost lattice site in
the original 1D real space lattice), we have $\varepsilon_{N^{2}}=0$,
indicating no local dissipation. In the Fock state lattice, this corresponds
to the site in the lower left corner (white dot). Along the edges originating
from the $l=N^{2}$-site, there is moderate dissipation (gray dots), while in
the bulk, the dissipation is most pronounced (black dots). Let us introduce
$\varepsilon_{b}$ and $\varepsilon_{e}$ for the bulk and edge on-site
“energies” respectively (the black and gray sites in the figure), and $t_{0}$
and $t_{d}$ for the horizontal/vertical (nearest neighbor) and diagonal (next
nearest neighbor) tunneling amplitudes respectively. These lattice parameters
are related to $\gamma$ as shown in Tab. 1. If there were no tunnelings
($t_{0}=t_{d}=0$), the system (ground/steady state) would localize at the site
represented by the white dot. The tunneling terms tend to delocalize the
steady state. As $\gamma$ increases, the dissipation-induced localization
becomes stronger, but at the same time, tunneling-induced delocalization also
strengthens. Before the critical point, $\gamma<\gamma_{c}=1/2$, the nearest
neighbor tunneling rates ($2|t_{0}|$) are greater than $\varepsilon_{b}$,
causing the system to delocalize. Beyond the critical point,
$\gamma>\gamma_{c}$, where $2|t_{0}|<\varepsilon_{b}$, the edge skin mode
wins. This reasoning demonstrates how the transition can be understood from
the Liouvillian Fock state lattice.

Nearest neighbour tunneling | $\hskip 5.69046ptt_{0}=\pm i(1-\gamma)\hskip 5.69046pt$
---|---
Diagonal tunneling | $t_{d}=2\gamma$
Bulk “energy” | $\varepsilon_{b}=-2\gamma$
Edge “energy” | $\varepsilon_{e}=-\gamma$
Table 1: Table of the parameters forming the Fock state lattice of Fig. 4,
i.e. green arrows – $t_{d}$, red arrows – $t_{0}$, black dots (bulk sites) –
$\varepsilon_{b}$, and gray dots (edge sites) – $\varepsilon_{e}$. Not shown
in the figure are the upper and right edges, where onsite “energies” are set
to $-3\gamma$. For sufficiently large lattices, these edge sites should not
significantly impact the criticality.

Figure 4: The Liouvillian Fock state lattice. Each lattice point represents a
Liouvillian Fock state $|l\rangle\\!\rangle\leftrightarrow|n\rangle\langle
m|$, and the colors denote the onsite (imaginary) “energy” with white ($=0$)
$>$ gray $>$ black. Consequently, apart from the lower left corner site, the
“population” within the other sites decays are suppressed by dissipation. The
arrows describe the non-vanishing tunneling elements. If quantum jumps, i.e.,
fluctuations, are disregarded, there would be no diagonal tunneling (green
arrows between next-nearest neighbors). The red arrows alone do not have a
preferred direction. The onsite dissipation tend to localize the state to the
lower left corner, while the tunneling terms counteract the localization. For
the HN model, one finds the same square lattice, but there is no diagonal
tunneling, and the horizontal/vertical tunnelings are imbalanced with
amplitudes $1\pm\gamma$, as expected for the HN model.

It is worth noting that similar delocalization-localization transitions have
been discussed in the past with a model described as

$\frac{\partial\hat{\rho}}{\partial
t}=i\left[\hat{\rho},\hat{S}_{x}\right]+\frac{\gamma}{S}\left(2\hat{S}^{-}\hat{\rho}\hat{S}^{+}-\hat{S}^{+}\hat{S}^{-}\hat{\rho}-\hat{\rho}\hat{S}^{+}\hat{S}^{-}\right).$
(38)

This model was introduced in the late 1970s to study quantum optical
bistability optbis . In this context, the $\hat{S}$-operators represent
collective spin operators, and $S$ is the total (conserved) spin. The
Hamiltonian can be seen as describing a classical drive of the spin, while the
Lindblad dissipation represents spontaneous decay of the spin. One attractive
aspect of this model is that the steady state can be determined analytically
optbis2 , and it exhibits a mean-field critical point at $\gamma_{c}=1/2$,
where the system transitions from being magnetized (localized) to paramagnetic
(delocalized) optbiscrit2 . The transition is continuous, without any apparent
spontaneous symmetry breaking, which has generated some debate optbiscrit2 .

To draw a connection to the model studied in this paper, it is important to
note that the Hamiltonian component $\hat{S}_{x}$ tends to delocalize the
state in the spin Fock basis (the $|S,m\rangle$-eigenstates of $\hat{S}_{z}$),
while the dissipative part drives the state toward the “edge” $|S,-m\rangle$.
We have numerically solved the LME using $\hat{H}=\hat{E}_{x}$ and
$\hat{L}=\hat{E}$, with $\hat{E}_{x}$ and $\hat{E}$ being operators from the
Euclidean algebra (21). This resulted in a similar phase transition as
observed for the LME (30). Consequently, we draw a comparison between two
models: one supporting an $SU(2)$ algebra and the other an Euclidean algebra,
both being otherwise equivalent. Of course, the spin operators come with a
square-root normalization factor when acting on the Fock states, but this
factor does not alter the Fock state lattice geometry; it induces a strain in
the lattice pil . It is important to mention that this connection is relevant
only for a finite lattice, as in the infinite case, all operators would
mutually commute, i.e.
$\left[\hat{E},\hat{E}^{\dagger}\right]=\left[\hat{E},\hat{E}_{x}\right]=\left[\hat{E}^{\dagger},\hat{E}_{x}\right]=0$.

The critical behavior of the bistability model has been extensively studied in
Ref. optbiscrit2 . It was argued that the critical behavior, devoid of
symmetry breaking, arises from the transition ‘softening’ of a first-order
transition into a continuous one. In the current model, we appear to observe
similar universal behavior. Specifically, the transition is continuous, the
steady state remains unique throughout, and hence, there is no apparent
symmetry breaking occurring. It is worth noting that quenched disorder can
alter the nature of a transition from first to second order in classical
critical models qd . Tri-critical points provide another example where the
order of a transition changes, and in the Potts model, the transition can
shift from first to second order as a system parameter varies potts .
Continuous phase transitions without symmetry breaking can also occur in
fermionic models when the Fermi sea undergoes volume changes fermisea .

It remains unclear whether the mechanism underlying the observed criticality
in this open system (and in the bistability model of Eq. (38)) differs in
nature from those listed in the references mentioned above. After all, we are
dealing with an open, non-Hamiltonian system. In Ref. optbiscrit2 , it was
noted that the full Hamiltonian model, including the degrees of freedom of the
environment, exhibited a first-order phase transition, and it was only in the
limiting case of infinite separation in time-scales that the initially first-
order, discontinuous phase transition became continuous. Something similar
might occur in our system, such as starting from the Hamiltonian (31) and
coupling it to a bath of oscillators.

###  IV.2 Non-hermitian sensors

In recent years, there have been several proposals on how systems described by
NH “Hamiltonians” can be leveraged to enhance sensor performance. Various
concepts for these implementations have been explored, including harnessing
the non-analyticity associated with exceptional points sensor3 , non-unitary
evolution sensor2 , and the phenomenon of exponential sensitivity sensor .
Motivated in part by a recent experimental demonstration sensor4 , our focus
will be on the latter aspect and its application to the HN model.

We aim to address two key questions in this context: firstly, what role do
fluctuations play in the sensor setup, and secondly, how does disorder impact
the performance of the NH sensor?

Figure 5: Logarithm of the HN auto-correlation Function (41), shown as a
function of the system size $N$ for four different asymmetry parameters
$\delta$ (solid black lines). In all four examples, the perturbation strength
is $\epsilon=10^{-10}$. The linearity of the curves clearly demonstrates the
exponential $N$-dependence and the point at which the perturbation $\delta\nu$
transitions from being real to complex, indicated as a “knee” in the curves.
The dashed red line represents an example of the disordered HN sensor
(averaged over 1000 disorder realizations) with a disorder strength $W=0.0005$
(i.e., $0.05\%$ of the tunneling rate). In the case of disorder, the sensor
must be of sufficient size for the perturbation to outweigh the effects of
random disorder. In practice, there exists a window, as marked for the
$\delta=0.25$ curve in the plot, within which the sensor operates effectively
(see also Fig. 7). For a sensor that is too small, the perturbation is
overshadowed by disorder, and for a sensor that is too large, the perturbation
becomes complex.

####  IV.2.1 Non-hermitian sensor with fluctuations

Let us begin by summarizing the fundamental concept of the NH sensor sensor .
As previously mentioned, the skin effect is directly linked to exponential
sensitivity. To illustrate this, let us consider the $N$-site open BC HN
model. The spectrum for this model was given in Eq. (22), and we observed how
it is real and symmetric around zero. This symmetry implies that for an odd
number of sites, there is a zero eigenvalue, denoted as
$\nu_{z}\equiv\nu_{j}=0$ for $j=(N+1)/2$. Now, let us further assume a very
weak (real) coupling denoted as $|\epsilon|\ll 1$, which connects the first
and last sites. This coupling represents a non-local perturbation of the form:

$\hat{V}=\epsilon\left(\hat{a}_{1}^{\dagger}\hat{a}_{N}+h.c.\right).$ (39)

Our goal here is to determine the value of $|\epsilon|$. When treating
$\hat{V}$ as a perturbation, we find (for NH, perturbation theory involves
both left and right unperturbed states since they form the orthogonal states
used for the resolution of identity) that the lowest-order corrections to the
zeroth eigenvalue scale as sensor

$\nu_{z}\rightarrow\nu_{z}=\delta\nu,\hskip 19.91684pt\delta\nu\sim\epsilon
e^{\alpha N}.$ (40)

Here, $\alpha$ depends on specific system details. Consequently, the
eigenvalue remains real but shifts away from zero. Notably, this shift can be
significant as long as $N$ is sufficiently large. Beyond a critical
perturbation $\epsilon_{c}(N)$ (which depends on the system size and other
system parameters), the scaling breaks down, and the eigenvalue becomes
complex. In line with the experimental work of Ref. sensor4 , if the system is
initially prepared in the zero eigenvalue state $|\varphi_{z}^{R}\rangle$ and
evolves for a brief time $t$ under the perturbed Hamiltonian
$\hat{H}_{\mathrm{HN}}+\hat{V}$, resulting in the state $|\psi\rangle$, then
the decay of the auto-correlation function

$A(N)=|\langle\varphi_{z}^{L}|\psi\rangle|$ (41)

serves as a measure of the perturbation.

Consequently, for a given perturbation strength $\epsilon$, due to the
exponential sensitivity, the logarithm of $A(N)$ should exhibit a linear
relationship with the system size N as long as $\epsilon<\epsilon_{c}(N)$. The
results, obtained from numerical time-propagation of the initial state, are
presented in Fig. 5 as solid black curves. In the figure, we display the
logarithm of the auto-correlation function for four different HN parameters,
denoted as $\delta$. It is worth recalling that in the Lindblad realization of
the HN model, $\delta$ corresponds to the rate $\gamma$. As expected, the
anticipated linear dependence on $N$ is demonstrated in all four examples
within the figure.

It is particularly noteworthy that for large $\delta$ values, the exponential
dependence as per Eq. (40) persists only for relatively small system sizes.
Conversely, this regime can extend over significantly larger system sizes for
small $\delta$ values. The slope of the curves is determined by the parameter
$\alpha$ in Eq. (40), which is found to scale as $\alpha\sim\delta$. This
relationship explains why rather large system sizes $N$ are required to reach
the critical $\epsilon_{c}$ for small $\delta$ values. In practical
implementations, state tomography provides the time-evolved state
$|\psi(t)\rangle$. With prior knowledge of $|\varphi_{z}^{L}\rangle$, it
becomes possible to estimate the decay of $A(N)$ and, consequently,
$\delta\nu$. This approach has been realized in classical light-pulses in
waveguides, where fluctuations likely play a less significant role (see
discussions in Sec. V) sensor4 .

Figure 6: The spectrum of the unperturbed (black dots) and perturbed
Liouvillian (red dots) for the localized phase (a) and the delocalized phase
(b). In the delocalized phase of (a) with $\delta=0.25$, it is challenging to
distinguish between red and black dots. However, in the localized phase (b)
with $\delta=0.75$, there is a clear exponential sensitivity. These results
are based on a system size of $N=61$ and a perturbation strength of
$\epsilon=10^{-10}$.

Now, turning our attention to the LME, we consider the same type of
perturbation (39). The exponential sensitivity in this context is a
consequence of the skin effect. As we have discussed, the skin effect can be
lost in the LME, leading to the steady state becoming delocalized when
$\gamma<\gamma_{c}=1/2$. Consequently, we expect that exponential sensitivity
is only present in the localized phase. Indeed, within this phase, the
eigenvectors of the Liouvillian become localized to the corner of the
Liouvillian Fock state lattice com2 . In Fig. 6, we provide two examples of
the Liouvillian spectrum: one in the delocalized phase (a) and the other in
the localized phase (b). As anticipated, the exponential sensitivity is no
longer present when the system transitions into the delocalized phase.

However, in the localized phase, the Liouvillian also exhibits exponential
sensitivity, a phenomenon that has been observed in the past nhlindblad1 ;
nhlindblad2 ; nhlindblad2b ; sofia . Nevertheless, the challenge remains to
identify an observable quantity that is exponentially sensitive and capable of
extracting information about the perturbation.

For the HN model, a natural choice is to consider the state associated with
the zero eigenvalue, $|\varphi_{z}^{R}\rangle\\!\rangle$. As we have seen,
this is due to the initial decay of its auto-correlation function, which
directly reflects the magnitude of the perturbation. In the case of the full
LME, the primary option would be the steady state, denoted as
$\hat{\rho}_{\mathrm{ss}}^{(\epsilon=0)}$ (where $\epsilon=0$ signifies the
unperturbed state). We can analyze how this steady state evolves under the
perturbed Liouvillian. Alternatively, we could also contemplate using the same
initial state $|\varphi_{z}^{R}\rangle$ as in the HN model and explore the
influence of fluctuations on the auto-correlation function (41). Through
numerical simulations, we find no signs of exponential sensitivity in either
of these scenarios.

In the latter case, the decay of the auto-correlation function is
predominantly driven by the relaxation towards the system’s steady state, and
it occurs on a timescale proportional to the inverse Liouvillian gap (12). The
perturbation-induced decay happens on an entirely different timescale and gets
obscured by the relaxation of the steady state. Extending the steady state
relaxation time by moving closer to the critical point does not result in a
favorable situation, as it rapidly suppresses exponential sensitivity when the
system becomes more delocalized. In our numerical experiments, we have not
identified a regime where the system effectively functions as a sensor when
initialized with the $|\varphi_{z}^{R}\rangle$ state.

This same argument applies to initializing the system in the unperturbed
steady state, meaning that the relaxation of
$\hat{\rho}_{\mathrm{ss}}^{(\epsilon=0)}$ to the final perturbed steady state
$\hat{\rho}_{\mathrm{ss}}^{(\epsilon)}$ dominates the evolution. Importantly,
this behavior does not exhibit a strong $N$-dependence. More precisely, we
find that the auto-correlation function
$A(N)=\mathrm{Tr}\left[\hat{\rho}_{\mathrm{ss}}^{(0)}\hat{\rho}_{\mathrm{ss}}^{(\epsilon)}\right]$
follows a linear $N$-dependence rather than an exponential one.

In summary, even though the spectrum of the Liouvillian exhibits exponential
sensitivity to non-local perturbations in the localized phase, it is not
evident how this can be effectively harnessed for sensing purposes.
Generalizing the concept of the HN sensor does not appear to yield favorable
results. It remains uncertain whether there might be another experimentally
measurable quantity that could salvage the sensor’s performance when
accounting for fluctuations.

####  IV.2.2 Non-hermitian sensor with disorder

Another significant limitation of the NH sensor that we need to address, which
does not stem from environment-induced fluctuations, pertains to disorder.
Consider the presence of local quenched disorder, which can arise in an
imperfect sensor, causing the actual Hamiltonian to take the form:

$\hat{H}_{\mathrm{dHN}}=\hat{H}_{\mathrm{HN}}+\sum_{n=1}^{N}\kappa_{n}\hat{n}_{n},$
(42)

where $\kappa_{n}\in\left[-W,W\right]$ represents a random onsite offset.
Here, $W$ denotes the disorder strength. We assume that $W\ll 1$, which means
it is orders of magnitude smaller than the tunneling rate, ensuring that the
system is not localized on any relevant length scales. However, it is
important to note that $W$ can be larger than the perturbation strength that
the sensor is meant to measure. Therefore, a potential breakdown of the sensor
should not be attributed to hindrance in propagation due to localization.

Figure 7: Gray regions indicate the parameter regimes in the $NW$-plane (note
the logarithmic scale of the disorder strength $W$) where the sensor can
effectively operate. Above this region, the $N$-dependence of the perturbation
is no longer exponential, and below, disorder-induced noise dominates the
detector signal. The perturbation strength is kept the same as in Figure 5,
i.e., $\epsilon=10^{-10}$.

In Fig. 5, we present the sensor’s performance in the presence of quenched
disorder, indicated by the red dashed line. The disorder strength is
approximately $0.05\%$ of the tunneling strength, and to reduce the scatter,
we averaged over 1000 disorder realizations. Our findings reveal that below a
certain system size, denoted as $N_{l}$, the sensor primarily detects the
disorder, with the perturbation signal getting lost in the noise generated by
the disorder. However, beyond $N_{l}$, the perturbation signal begins to
dominate the disorder noise, thanks to its exponential increase, and the
disorder’s impact on the signal diminishes significantly.

It is important to note that above a certain upper system size, $N_{u}$, the
perturbation signal no longer follows the exponential form described in Eq.
(40). Consequently, when given a perturbation $\epsilon$ and a disorder
strength $W$, there exists a window of system sizes $N_{l}<N<N_{u}$ within
which the NH sensor operates effectively. These windows, corresponding to the
same values of $\delta$ as seen in Fig. 5, are displayed in Fig. 7.

In summary, when $\delta$ is small, it implies the need for larger sensor
sizes $N$. In practice, achieving extremely large chains may not be feasible,
making it desirable to use a larger $\delta$. Our results are computed using
the same methodology as for the non-disordered sensor, involving the
propagation of the initial state and the analysis of the auto-correlation
function’s decay. We have also verified numerically that similar results can
be obtained by directly extracting $\delta\nu$ from the full spectrum, rather
than indirectly calculating it from $A(N)$.

##  V Discussion and Concluding remarks

While NH Hamiltonians offer utility in modeling a quantum system’s interaction
with its environment, they can also lead to complex, non-physical outcomes. A
significant challenge in this endeavor arises from the intricate entanglement
between the quantum system and its environment, resulting in a mixed state
representation rather than the pure state representation
$\hat{\rho}=|\psi\rangle\langle\psi|$. However, it is worth noting that there
are scenarios in which the quantum system’s state remains nearly pure. In such
cases, relevant observables
$\mathcal{O}=\mathrm{Tr}\left[\hat{\rho}\hat{O}\right]$ can be well described
by a pure state. Many experimental activities involve classical emulations of
quantum systems, where classical systems obey equations of motion similar, or
equivalent, to those found in quantum systems described by NH Hamiltonians.
While classical systems are not represented by states in a Hilbert space, the
fluctuation theorem extends to both classical and quantum systems. In a strict
sense, any open system is subject to fluctuations from its environment, with
Brownian motion serving as a classic example. However, this influence may be
negligible for macroscopic objects. Classical states, described as coherent
states in the quantum realm, are known to be robust against fluctuations, such
as the state of the electromagnetic field originating from a laser com3 .

The role of fluctuations becomes a more intricate matter when dealing with
systems deep in the quantum regime. To circumvent this issue, researchers have
employed the concept of post-selection postselect , where experiments are
conducted under full observation, and only the data from experimental runs
that do not undergo “quantum jumps”, such as spontaneous photon emission, are
considered unravel . However, this approach often leads to a significantly
reduced probability of successful experimental runs over time, as they become
exponentially suppressed.

In this study, we explored the time evolution generated by the Liouvillian
without invoking measurement-induced projections or assuming a semi-classical
regime. Specifically, we focused on two well-studied phenomena within the
framework of NH QM: the skin effect and NH sensors.

While previous arguments suggested that both the skin effect and NH sensors
should persist when fluctuations are taken into account nhlindblad1 ;
nhlindblad2 ; nhlindblad4 , our study unveiled subtleties in this context,
leading to new insights. Importantly, the Liouvillian is not uniquely
determined by a NH “Hamiltonian”. This multiplicity is akin to the
purification of mixed states nc , where infinitely many different pure states
can construct a given mixed state. Consequently, there exist infinitely many
Liouvillians corresponding to equivalent NH “Hamiltonians”, and these
Liouvillians can exhibit qualitative differences.

In our study, we chose to investigate a Liouvillian represented by collective
quantum jumps, as defined in Eq. (30). We also considered an alternative model
with local jump operators, as defined in Eq. (27). Upon numerical analysis of
the latter model, we found that it lacks critical behavior – the delocalized
phase does not emerge in this case. This finding implies that the
collectiveness of quantum jumps is essential for the appearance of a critical
point. The criticality observed in our model is an example of a driven-
dissipative non-equilibrium phase transition lcrit . This specific type of
transition is qualitatively different from phase transitions described by the
Ginzburg-Landau paradigm cardy . The transition is continuous and lacks any
apparent symmetry breaking, yet the critical exponent for the correlation
length aligns with the characteristics of a mean-field transition. This
observation resonates with related models, such as the one represented in Eq.
(38), which has been explored in the context of optical bistability.

The existence of a delocalized phase signifies the breakdown of the NH skin
effect. In this phase, the system exhibits an extended steady state that
approximately populates the lattice sites uniformly. The coherences between
the sites practically vanish, implying that the steady state approximates a
maximally mixed state or an infinite-temperature state. In the case of
periodic BC, the maximally mixed state is the exact steady state, a concept
that was briefly mentioned in a previous study hatanonelsonlindblad , which
investigated a quadratic Liouvillian.

It is essential to recognize that not only the steady state becomes
delocalized for $\gamma<\gamma_{c}$, but the other Liouvillian eigenvectors
$\hat{\rho}_{j}$ also undergo delocalization in this phase. This
delocalization affects the entire spectrum of the Liouvillian.

Indeed, the absence of a skin effect in the delocalized phase results in the
system’s lack of exponential sensitivity. Consequently, it cannot be employed
for sensing purposes in this regime. However, when the losses are substantial,
i.e., $\gamma>\gamma_{c}$, the Liouvillian does exhibit exponential
sensitivity once it enters the localized phase.

Nonetheless, despite the presence of exponential sensitivity in this regime,
our exploration did not yield a suitable measure that could effectively
extract the quantity to be detected. Attempts to directly generalize
experiments like the one reported in sensor4 , which used classical light,
were ineffective in the context of Liouvillian systems. This inefficacy is
primarily due to the dominance of other mechanisms, such as relaxation towards
a steady state, during the early stages of evolution. We also considered
observing the long-term evolution and the fidelities of the resulting steady
states, but these did not exhibit exponential sensitivity either.

While the specific observables explored in this study did not yield the
desired results, it remains a possibility that more refined observables could
be identified for use in Liouvillian-based sensors. However, further
investigation and research would be required to identify and assess these
potential observables effectively.

In our study, we have also demonstrated that disorder introduces limitations
to the applicability of NH sensors, even without considering fluctuations.
This finding underscores the importance of the sensor’s size relative to the
disorder strength. The sensor must be sufficiently large to overcome the noise
due to disorder.

In particular, we observed that for the NH sensor to effectively detect a
perturbation signal, the signal strength must exceed a critical value that is
determined by the disorder strength. This means that the sensor must be
sufficiently large to overcome the detector noise generated by the disorder.

These findings emphasize the importance of carefully considering disorder
effects and ensuring that the sensor size and sensitivity are suitable for
practical applications. It also suggests that engineering sensors with greater
robustness against disorder may be necessary for reliable measurements in
realistic experimental conditions.

While our analysis has been centered on the HN model, the applicability of our
results to other models remains an open question. In Appendix B, however, we
extend our findings to the NH SSH chain and demonstrate that the qualitative
conclusions applies. Nevertheless, one can imagine the potential for further
research in other directions that addresses slightly different questions,
including topology, localization, and the realm of many-body NH physics. We
propose that the introduction of the Liouvillian Fock state lattice, as
illustrated in Fig. 4, offers a promising tool for gaining fresh insights.
Particularly fascinating is the observation that, for the 1D HN model, the
dimension of the Liouvillian Fock state lattice stands at 2D. This observation
leads to the intriguing speculation that, for a D-dimensional model coupled to
an environment, the relevant physics might manifest in D+1 dimensions.
Furthermore, the interplay of symmetries in NH quantum mechanics compared to
the full LME lsym  remains an open question, partially due to the complexity
of relating biorthogonal QM to real physical systems. These questions are left
unexplored for future studies.

##  Appendix A Open quantum systems

Figure 8: Description of the open system evolution. The (sub)system
$\mathcal{S}$ is represented by a state $\hat{\rho}_{\mathrm{s}}(t)$,
schematically depicted in blue in (a). The system interacts weakly with an
environment $\mathcal{E}$ (light blue), allowing information within
$\mathcal{S}$ to irreversibly dissipate into $\mathcal{E}$ with a
characteristic decay rate $\gamma$. In (b), the concept of the CPTP dynamical
map is presented. The full system comprises $\mathcal{S}$ and $\mathcal{E}$,
and its state $\hat{\rho}(t)$ undergoes unitary evolution under a Hamiltonian
$\hat{H}$, as indicated by the upper horizontal double arrow. The system state
$\hat{\rho}_{\mathrm{s}}(t)$ is derived by taking the partial trace of the
full state over the environment degrees of freedom, represented by vertical
downward arrows. The dynamical map $\hat{\rho}_{\mathrm{s}}(t)$, shown as the
lower horizontal arrow, governs the evolution of the system state,
transforming $\hat{\rho}_{\mathrm{s}}(0)$ into $\hat{\rho}_{\mathrm{s}}(t)$.
While the full state evolution is unitary and thus reversible, all other steps
(partial traces and dynamical map) are irreversible.

In this appendix, we outline some general concepts related to the dynamics of
open quantum systems. The typical scenario is illustrated in Fig. 1. In (a),
we depict a small system $\mathcal{S}$ interacting with a larger environment
$\mathcal{E}$. Information is exchanged between these subsystems, with the
rate of exchange determined by the parameter $\gamma$. This exchange includes
processes such as particle or energy losses and decoherence. The combined
system’s evolution is described by the Schrödinger equation. For instance, an
initial pure separable state evolves as

$|\Psi(t)\rangle=\hat{U}(t)|\psi_{\mathrm{s}}(0)\rangle\otimes|\psi_{\mathrm{E}}(0)\rangle.$
(43)

Here, the time-evolution operator $\hat{U}(t)$ is generated from the full
Hamiltonian
$\hat{H}=\hat{H}_{\mathrm{s}}+\hat{H}_{\mathrm{E}}+\hat{H}_{\mathrm{sE}}$,
where the first two terms represent the system and environment sub-
Hamiltonians, and the last term accounts for their interaction. The state of
the system at time $t$ is obtained by taking the partial trace of the full
state, as illustrated in Fig. 9 (b).

As time progresses, the initially separable state of the system-environment
becomes entangled. Given that the initial state is pure, this entanglement
between the system and the environment is reflected in the system’s reduced
state, $\hat{\rho}_{\mathrm{s}}(t)$, being mixed. For instance, the von
Neumann entropy
$S_{\mathrm{vN}}=-\mathrm{Tr}\left[\hat{\rho}(t)\ln\hat{\rho}(t)\right]>0$ is
used as a measure of the amount of entanglement [33]. From this point forward,
we will omit the subscript “s” for the reduced density operator of the system.

In the Markovian approximation, as assumed for the LME (2), this information
flow out of the system is forever lost to the environment bp . Furthermore, in
deriving equation (2), we also assume the validity of the Born approximation,
which implies that, due to the substantial difference in system sizes, the
state of the environment is unaffected by the presence of the small system bp
.

As previously mentioned, the loss of information can happen through either
dissipation or decoherence. In the former, we often think of particle losses,
while in the latter, it typically results from uncontrolled energy shifts.
Regardless of the specific mechanism, both processes tend to cause the state
of the system to become mixed in most cases. This concept is at the core of
fluctuation-dissipation or quantum regression theorems bp ; mw ; carmichael .
However, if one has complete access to the environment, in principle, it is
possible to extract all the information about the system. In this scenario,
the system’s state can be ascribed a pure state. This would be the case of
post-selection postselect , briefly touched upon in the final Sec. V.

Figure 9: Demonstration of the localized-delocalized transition for the NH SSH
model. In (a), we present the same information as Fig. 2 (a) for the HN model.
For large $\gamma$ values, the steady state becomes localized at the edge of
the chain, whereas for weak $\gamma$ values, it becomes delocalized over the
entire lattice, as evident in (b), which displays the width (45). The dashed
red curves correspond to $\chi=0.5$, and the solid black curves to $\chi=0.1$.
The dotted line in (b) represents the width for a maximally delocalized state.
It is evident that the SSH parameter $\chi$ smoothens the transitions, but
criticality appears to persist. This is further confirmed in (c) and (d) for
the real parts of the Liouvillian spectra.

##  Appendix B Liouvillian for non-hermitian SSH model

It is straightforard to generalize the HN to the Su–Schrieffer–Heeger (SSH)
model ssh . The Hamiltonian of the hermitian SSH model is

$\hat{H}_{\mathrm{SSH}}=(1-\chi)\sum_{n=1}^{N}\hat{a}_{n}^{\dagger}\hat{b}_{n}+(1+\chi)\sum_{n=1}^{N}\hat{a}_{n+1}^{\dagger}\hat{b}_{n}+h.c.$
(44)

The parameter $\chi$, which takes values in the interval $[-1,0]$, dictates
the relative strengths of consecutive tunneling amplitudes in the SSH model.
For $\chi=0$, the model reverts to the conventional tight-binding model.
However, for $\chi\neq 0$, the lattice exhibits a bipartite structure, with
each unit cell containing two distinct types of sites: the $a$-sites and the
$b$-sites. In cases with an odd number of sites, the SSH model supports a
zero-energy topological edge state, also known as symmetry-protected state,
which exhibit exponential localization near the chain’s edge. The specific
value of $\chi$ determines the size of the energy gap that separates these
edge states from the remaining bulk states, as well as the amount of
localization of the edge state. These edge states tend to be more isolated
from the rest, potentially enhancing their robustness for various
applications. Consequently, it is worthwhile to investigate how the findings
of the present paper extend to the SSH model. Transitioning to the SSH model
from the tight-binding Hamiltonian in Eq. (30) is straightforward. The results
of numerical simulations are illustrated in Fig. 9, displaying the spectra,
scaled positions (32), and widths

$\Delta\xi=\sqrt{\langle\hat{S}^{2}\rangle-\langle\hat{S}\rangle^{2}}.$ (45)

The plot clearly demonstrates that $\chi$ influences the transitions, shifting
them toward smaller values of $\gamma$, and the transition itself becomes
smoother. Nonetheless, based on the available numerical data, it appears that
the transition remains distinct and is not merely a crossover. A noteworthy
observation is that, even when the steady state is not precisely centered in
the middle of the lattice, its width remains close to the maximum.
Importantly, it should be noted that the model exhibits asymmetry concerning
the sign of $\chi$. This asymmetry arises from the drift induced in the
lattice by the jump operators, and reversing the sign of $\chi$ results in a
distinct model.

###### Acknowledgements.

The author acknowledges financial support from VR-Vetenskapsråset (The Swedish
Research Council), and is thankful for fruitful discussions with Elisabet
Edvardsson and Emil Bergholtz.

## References

  * (1) M. A. Miri, and A. Alu, Science, 363, eaar7709 (2019); Y. Ashida, Z. Gong, and M. Ueda, Adv. Phys. 69, 249 (2020); E. J. Bergholtz, J. C. Budich, and F. K. Kunst, Rev. Mod. Phys. 93, 015005 (2021).
  * (2) C. M. Bender and S. Boettcher, Phys. Rev. Lett. 80, 5243 (1998); C. M. Bender, Rep. Prog. Phys. 70, 947 (2007).
  * (3) C. M. Bender, D. C. Brody, and H. F. Jones, Phys. Rev. Lett. 89, 270401 (2002); D. C. Brody, J. Phys. A Math. Theo. 47, 035305 (2013); F. K. Kunst, Kunst, F. K., E. Edvardsson, J. C. Budich, and E. J. Bergholtz, Phys. Rev. Lett. 121, 026808 (2018).
  * (4) C. M. Bender, D. C Brody, and H. F Jones, Am. J. Phys. 11, 1095 (2003).
  * (5) S. Yao and Z. Wang, Phys. Rev. Lett. 121, 086803 (2018).
  * (6) N. Okuma, K. Kawabata, K. Shiozaki, and M. Sato, Phys. Rev. Lett. 124, 086801 (2020); X. Zhang, T. Zhang, M. H. Lu, and Y. F. Chen, Adv. Phys. X 7, 2109431 (2022).
  * (7) K. Kawabata, T. Numasawa, and S. Ryu, Phys. Rev. X 13, 021007 (2022).
  * (8) E. Edvardsson and E. Ardonne, Phys. Rev. B 106, 115107 (2022).
  * (9) R. Koch, and J. C. Budich, Euro. Phys. J. D 74, 1 (2020); C.-X. Guo, C.-H. Liu, X.-M. Zhao, Y. Liu, and S. Chen, Phys. Rev. Lett. 127, 116801 (2021).
  * (10) J. C. Budich and E. J. Bergholtz, Phys. Rev. Lett. 125, 180403 (2020).
  * (11) M. Parto, C. Leefmans, J. Williams, and A. Marandi, arXiv:2305.03282.
  * (12) R. Kubo, Rep. Prog. Phys. 29, 255 (1966); L. Mandel, and E. Wolf, Optical coherence and quantum optics, (Cambridge university press, Cambridge, 1995).
  * (13) H. J. Carmichael, Statistical Methods in Quantum Optics 1, (Springer, 2002).
  * (14) Y.-C. Lee, M.-H. Hsieh, S. T. Flammia, and R.-K. Lee, Phys. Rev. Lett. 112, 130404 (2014); J.-S. TANG, et al., Nature Phot. 10, 642 (2016); A. Kumari and U. Sen, arXiv:2202.02744.
  * (15) X. Zhan, K. Wang, L. Xiao, Z. Bian, Y. Zhang, B. C. Sanders, C. Zhang, and P. Xue, Phys. Rev. A 101, 010302(R) (2020).
  * (16) Y. Ashida and M. Ueda, Phys. Rev. Lett. 120, 185301 (2018); N. Matsumoto, K. Kawabata, Y. Ashida, S. Furukawa, and M. Ueda, Phys. Rev. Lett. 125, 260601 (2020); B. Barch, N. Anand, J. Marshall, E. Rieffel, and P. Zanardi, arXiv:2305.12054.
  * (17) S.-L. Chen, G.-Y. Chen, and Y.-N. Chen, Phys. Rev. A 90, 054301 (2014); A. K. Pati, arXiv:1404.6166; J.-W. Wen, C. Zheng, X.-Y. Kong, S.-J. Wei, T. Xin, and G.-L. Long, Phys. Rev. A 99, 062122 (2019).
  * (18) Y. T. Tu, Y. C. Tzeng, and P. Y. Chang, SciPost Phys. 12, 194 (2022); M. Fossati, F. Ares, and P. Calabrese, Phys. Rev. B 107, 205153 (2023); C. T. Hsieh, and P. Y. Chang, 6, 062 (2023)
  * (19) N. Hatano and N. D. Nelson, Phys. Rev. Lett. 77, 570 (1996); ibid, Phys. Rev. B 56, 8651 (1997).
  * (20) R. Bonifacio and L. A. Lugiato, Phys. Rev. A 18, 1129 (1978); R. Bonifacio and L. A. Lugiato, Phys. Rev. Lett. 40, 1023 (1978); G. P. Agrawal and H. J. Carmichael, Phys. Rev. A 19, 2074 (1979); P. D. Drummond and D. F. Walls, J. Phys. A 13, 725 (1980); S. R. K. Rodriguez, W. Casteels, F. Storme, N. Carlon Zambon, I. Sagnes, L. Le Gratiet, E. Galopin, A. Lemaître, A. Amo, C. Ciuti, and J. Bloch, Phys. Rev. Lett. 118, 247402 (2017).
  * (21) R. R. Puri and S. V. Lawande, Phys. Lett. A 72, 200 (1979); S. V. Lawande, R. R. Puri, and S. S. Hassan, J. Phys. B 14, 4171 (1981).
  * (22) K. G. Zloshchastiev, and A. Sergi, J. Mod. Opt. 61, 1298 (2014).
  * (23) F. Minganti, A. Miranowicz, R. W. Chhajlany, and F. Nori, Phys. Rev. A, 100, 062131 (2019).
  * (24) A. McDonald, R. Hanai, and A. A. Clerk, Phys. Rev. B 105, 064302 (2022).
  * (25) S. Longhi, Phys. Rev. B 102, 201103(R) (2020).
  * (26) F. Song, S. Yao, and Z. Wang, Phys. Rev. Lett. 123, 170401 (2019); C.-H. Liu, K. Zhang, Z. Yang, and S. Chen, Phys. Rev. Research 2, 043167 (2020); X. Li, M. A. Begaowe, S. Zhang, and B. Flebus, arXiv:2307.15792.
  * (27) N. Okuma and M. Sato, Phys. Rev. B 103, 085428 (2021); F. Yang, Q.-D. Jiang, and E. J. Bergholtz, Phys. Rev. Research 4, 023160 (2022); Z. Zhou and Z. Yu, Phys. Rev. A 106, 032216 (2022).
  * (28) F. Koch and J. C. Budich, Phys. Rev. Research 4, 013113 (2022).
  * (29) A. McDonald and A. A. Clerk, Nature Comm. 11, 5382 (2020); L. Bao, B. Qi, D. Dong, and F. Nori, Phys. Rev. A 103, 042418 (2021); L. Bao, B. Qi, F. Nori, and D. Dong, arXiv:2303.16575.
  * (30) X. Niu, J. Li, S. L. Wu, and X. X. Yi, arXiv:2202.12591.
  * (31) M. A. Nielsen and I. L. Chuang, Quantum Computation and Quantum Information, (Cambridge University Press, 2000).
  * (32) J. Preskill, Lecture notes for physics 229: Quantum information and computation, California Institute of Technology 16, 1 (1998); V. Vedral, Introduction to Quantum Information Science, (Oxford University Press, 2006).
  * (33) I. Bengtsson, Göran Lindblad in memoriam, arXiv:2307.10621.
  * (34) G. Lindblad, Commun. Math. Phys. 48, 119 (1976).
  * (35) H.-P. Breuer and F. Petruccione, The Theory of Open Quantum Systems, (Oxford University Press, 2002).
  * (36) S. Clark, A. Peng, M. Gu, and S. Parkins, Phys. Rev. Lett. 91, 177901 (2003); J. T. Barreiro, P. Schindler, O. Gühne, T. Monz, M. Chwalla, C. F. Roos, M. Hennrich, and R. Blatt, Nature Phys. 6, 943 (2010); J. Cho, S. Bose, and M. Kim, Phys. Rev. Lett. 106, 020504 (2011); S. Diehl, E. Rico, M. A. Baranov, and P. Zoller, Nature Phys. 7, 971 (2011); M. Müller, S. Diehl, G. Pupillo, P. Zoller, et al., Adv. At. Mol. Opt. Phys. 61, 1 (2012); P. Schindler, M. Müller, D. Nigg, J. T. Barreiro, E. A. Martinez, M. Hennrich, T. Monz, S. Diehl, P. Zoller, and R. Blatt, Nature Phys. 9, 361 (2013).
  * (37) H. Spohn, Lett. Math. Phys. 2, 33 (1977); F. Alberto, Comm. Mth. Phys. 63, 269 (1978); D. Nigro, J. Stat. Mech. 4, 043202 (2019); H. Yoshida, arXiv:2309.00335.
  * (38) V. V. Albert and L. Jiang, Phys. Rev. A 89, 022118 (2014).
  * (39) P. Hedvall and J. Larson, arXiv: 1712.01560.
  * (40) J. A. Gyamfi, Eur. J. Phys. 41, 063002 (2020).
  * (41) J. Larson and S. Qvarfort, Open Sys. Quant. Inf. 30, 2350008 (2023).
  * (42) N. G. van Kampen, Stochastic processes in physics and chemistry, (Elsevier, 1992).
  * (43) P Alsing and H J Carmichael, Quant. Opt. 3 13 (1991); J. Dalibard, Yvan Castin, and Klaus Mølmer Phys. Rev. Lett. 68, 580 (1992); A. J. Daley, Adv. Phys. 63, 77 (2014).
  * (44) A. Rivas and S. F. Huelga, Open Quantum Systems, (Springer, 2012).
  * (45) E. M. Kessler, G. Giedke, A. Imamoglu, S. F. Yelin, M. D. Lukin, and J. I. Cirac, Phys.l Rev. A 86, 012116 (2012); F. Minganti, A. Biella, N. Bartolo, and C. Ciuti, Phys. Rev. A 98, 042118 (2018).
  * (46) E. M. Kessler, G. Giedke, A. Imamoglu, S. F. Yelin, M. D. Lukin, and J. I. Cirac, Phys. Rev. A 86, 012116 (2012); F. Minganti, A. Biella, N. Bartolo, and C. Ciuti, Phys. Rev. A 98, 042118 (2018).
  * (47) T. Mori and T. Shirai, Phys. Rev. Lett. 125, 230604 (2020); T. Haga, M. Nakagawa, R. Hamazaki, and M. Ueda, Phys. Rev. Lett. 127, 070402 (2021).
  * (48) W. Casteels, R. Fazio, and C. Ciuti, Phys. Rev. A 95, 012128 (2017); F. Vicentini, F. Minganti, R. Rota, G. Orso, and C. Ciuti, Phys. Rev. A 97, 013853 (2018); R. Rota, F. Minganti, C. Ciuti, and V. Savona, Phys. Rev. Lett. 122, 110405 (2019); M. Nakagawa, N. Kawakami, and M. Ueda, Phys. Rev. Lett. 126, 110404 (2021); G. T. Landi, D. Poletti, and G. Schaller, Rev. Mod. Phys. 94, 045006 (2022).
  * (49) B. Buca, J. Tindall, and D. Jaksch, Nature Commun. 10, 1730 (2019).
  * (50) W. D. Heiss, J. Phys. A Math. Theo. 45, 444016 (2012).
  * (51) S. S. Roy, S. Bandyopadhyay, R. C. de Almeida, and P. Hauke, arXiv:2309.00049.
  * (52) C. Y. Ju, A. Miranowicz, G. Y. Chen, and F. Nori, Phys. Rev. A 100, 062118 (2019); S. S. Bhosale, B. Rath, and P. K. Panigrahi, Quant. Rep. 3, 417 (2021); D. C. Brody, J. Phys. A: Math. Theor. 49, 10LT03 (2016).
  * (53) F. K. Kunst, E. Edvardsson, J. C. Budich, and E. J. Bergholtz, Phys. Rev. Lett. 121, 026808 (2018).
  * (54) F. G. Scholtz, H. B. Geyer, and F. J. W. Hahne, Ann. Phys. 213, 74 (1992); A. Mostafazadeh, J. Math. Phys. 44, 974 (2003); E. Edvardsson, J. L. K. König, and M. Stålhammar, arXiv:2212.06004.
  * (55) Motivating NH QM as emerging from a more complete description can be done in various ways, and starting from the LME is only one of them. This is, however, the one approach we discuss in this paper.
  * (56) J.M. Torres, Phys. Rev. A 89, 052133 (2014).
  * (57) T. Prosen, New J. Phys. 10, 043026 (2008).
  * (58) N. Hatano and D. R. Nelson, Phys. Rev. B 58, 8384 (1998).
  * (59) A. B. Klimov and S. M. Chumakov. A Group-theoretical Approach to Quantum Optics: Models of Atom-field Interactions, (John Wiley and Sons, 2009).
  * (60) F. Roccati, Phys. Rev. A 104, 022215 (2021).
  * (61) M. Lewenstein, A. Sanpera, V. Ahufinger, B. Damski, A. Sen, and U. Sen, Adv. Phys. 56, 243 (2007); I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80, 885 (2008).
  * (62) Z. Wang, Y. Lu, Y. Peng, R. Qi, Y. Wang, and J. Jie Phys. Rev. B 108, 054313 (2023).
  * (63) C. Gardiner, and P. Zoller, The quantum world of ultra-cold atoms and light book II: the physics of quantum-optical devices, (World Scientific Publishing, 2015); J. Larson and T. K. Mavrogordatos, The Jaynes-Cummings model and its descendants, (IOP Publishing, 2022).
  * (64) W. Chen, S. Kaya Özdemir, G. Zhao, J. Wiersig, and L. Yang, Nature 548, 192 (2017); H. Hodaei, A. U. Hassan, S. Wittek, H. Garcia-Gracia, R. El-Ganainy, D. N. Christodoulides, and M Khajavikhan. Nature 548, 187 (2017); H.-K. Lau and A. A. Clerk, Nature Comm. 9, 4320 (2018); M. Zhang, W. Sweeney, C. W. Hsu, L. Yang, A. D. Stone, and L. Jiang, Phys. Rev. Lett. 123, 180501 (2019);M. Parto, Y. G. Liu, B. Bahari, M. Khajavikhan, and D. N. Christodoulides, Nanophot. 10, 403 (2020); S. Soleymani, Q. Zhong, M. Mokim, S. Rotter, R. El-Ganainy, and S. K. Özdemir, Nature Comm. 13, 599 (2022).
  * (65) J. Cardy, Scaling and renormalization in statistical physics, (Cambridge University Press, 1996).
  * (66) P. Ivanov, N. Vitanov, and K. Bergmann, Phys. Re- v. A 70, 063409 (2004); M. S. Sarandy, D. A. Lidar, Phys. Rev. Lett. 95, 250503 (2005); T. Mathisen and J. Larson, Entropy 20, 20 (2018).
  * (67) H. Cai and D.-W. Wang, National Sci. Rev. 8, nwaa196 (2021); P. M. Saugmann and J. Larson, arXiv:2203.13813.
  * (68) D. F. Walls, P. D. Drummond, S. S. Hassan, and H. J. Carmichael, Prog. Theor. Phys. 64, 307 (1978); P. D. Drummond and H. J. Carmichael, Opt. Commun. 27, 160 (1978); P. D. Drummond, Phys. Rev. A 22, 1179 (1980); S. Schneider and G. J. Milburn, Phys. Rev. A 65, 042107 (2002); A. Dombi, A. Vukics, and P. Domokos, J. Phys. B 46, 224010 (2013); J. J. Mendoza-Arenas, S. R. Clark, S. Felicetti, G. Romero, E. Solano, D. G. Angelakis, and D. Jaksch, Phys. Rev. A 93, 023821 (2016); R. M. Wilson, K. W. Mahmud, A. Hu, A. V. Gorshkov, M. Hafezi, and M. Foss-Feig, Phys. Rev. A 94, 033801 (2016); W. Casteels, R. Fazio, and C. Ciuti, Phys. Rev. A 95, 012128 (2017); M. Foss-Feig, P. Niroula, J. T. Young, M. Hafezi, A. V. Gorshkov, R. M. Wilson, and M. F. Maghrebi, Phys. Rev. A 95, 043826 (2017)
  * (69) J. Hannukainen and J. Larson, Phys. Rev. A 98, 042113 (2018).
  * (70) M. Aizenman and J. Wehr, Phys. Rev. Lett. 62, 2503 (1989).
  * (71) B. Nienhuis, A. N. Berker, E. K. Riedel, and M. Schick, Phys. Rev. Lett. 43, 737 (1979); E. Domany, M. Schick, and R. H. Swendsen, Phys. Rev. Lett. 52, 1535 (1984).
  * (72) T. Senthil, Matthias Vojta, and Subir Sachdev Phys. Rev. B 69, 035111 (2004); Q. Si, and F. Steglich, Science, 329, 1161 (2010); N. Maksimovic et al., Science 375, 76 (2022).
  * (73) The degree of localization in the Liouvillian Fock state lattice can be calculated from the inverse partition ratio, e.g. $IPR_{j}=1/\sum_{l}|\langle\\!\langle l|\varphi_{j}^{R}\rangle\\!\rangle|^{4}$, where $|l\rangle\\!\rangle$ is a Liouvillian Fock state, and $|\varphi_{j}^{R}\rangle\\!\rangle$ is the normalized $j$’th eigenvector of the Liouvillian. One finds (not shown in this paper) that all eigenvectors are to some high localized for $\gamma>\gamma_{c}$, and they are much less localized for $\gamma<\gamma_{c}$ (even though the limiting value of $IPR_{j}=N^{2}$ is never reached even in the delocalized phase).
  * (74) Photon losses of the laser resonator can be mimicked by a Lindblad jump operator $\hat{L}=\hat{a}$, where $\hat{a}$ is the photon annihilation operator of the light mode. Since a coherent state is an eigenstate of this operator, its action on the state does not alter it; $\hat{a}|\alpha\rangle=\alpha|\alpha\rangle$. As soon as the state is not a coherent one, the presence of fluctuations will influence it.
  * (75) M. Naghiloo, M. Abbasi, Y. N. Joglekar, and K. W. Murch, Nature Physics, 15, 1232 (2019); G. L. Zhang, D. Liu, H. M. Yung, Scientific Rep. 11, 13795 (2021); A. Quinn, J. Metzner, J. E. Muldoon, I. D. Moore, S. Brudney, S. Das, D. T. C. Allcock and Y. N. Joglekar, arXiv:2304.12413; S. Erdamar, M. Abbasi B. Ha W. Chen, J. Muldoon Y. Joglekar and K. W. Murch, arXiv:2309.12393.
  * (76) P. Alsing, H. J. Carmichael, J. Euro. Opt. Soc. B 3, 13 (1991); J. Dalibard, Y. Castin, and K. Mölmer, Phys. Rev. Lett. 68, 580 (1992).
  * (77) V. V. Albert and L. Jiang, Phys. Rev. A 89, 022118 (2014).
  * (78) S. Lieu, Phys. Rev. B 97, 045106 (2018); S. Yao and Z. Wang, Phys. Rev. Lett. 121, 086803 (2018).