font=small

††thanks<EMAIL_ADDRESS>
# Reaction-limited quantum reaction-diffusion dynamics

Gabriele Perfetto  Institut für Theoretische Physik, Eberhard Karls
Universität Tübingen, Auf der Morgenstelle 14, 72076 Tübingen, Germany.
Federico Carollo  Institut für Theoretische Physik, Eberhard Karls Universität
Tübingen, Auf der Morgenstelle 14, 72076 Tübingen, Germany.  Juan P. Garrahan
School of Physics, Astronomy, University of Nottingham, Nottingham, NG7 2RD,
UK.  Centre for the Mathematics, Theoretical Physics of Quantum Non-
Equilibrium Systems, University of Nottingham, Nottingham, NG7 2RD, UK  Igor
Lesanovsky  Institut für Theoretische Physik, Eberhard Karls Universität
Tübingen, Auf der Morgenstelle 14, 72076 Tübingen, Germany.  School of
Physics, Astronomy, University of Nottingham, Nottingham, NG7 2RD, UK.  Centre
for the Mathematics, Theoretical Physics of Quantum Non-Equilibrium Systems,
University of Nottingham, Nottingham, NG7 2RD, UK

###### Abstract

We consider the quantum nonequilibrium dynamics of systems where fermionic
particles coherently hop on a one-dimensional lattice and are subject to
dissipative processes analogous to those of classical reaction-diffusion
models. Particles can either annihilate in pairs, $A+A\to\emptyset$, or
coagulate upon contact, $A+A\to A$, and possibly also branch, $A\to A+A$. In
classical settings, the interplay between these processes and particle
diffusion leads to critical dynamics as well as to absorbing-state phase
transitions. Here, we analyze the impact of coherent hopping and of quantum
superposition, focusing on the so-called reaction-limited regime. Here,
spatial density fluctuations are quickly smoothed out due to fast hopping,
which for classical systems is described by a mean-field approach. By
exploiting the time-dependent generalized Gibbs ensemble method, we
demonstrate that quantum coherence and destructive interference play a crucial
role in these systems and are responsible for the emergence of locally
protected dark states and collective behavior beyond mean-field. This can
manifest both at stationarity and during the relaxation dynamics. Our
analytical results highlight fundamental differences between classical
nonequilibrium dynamics and their quantum counterpart and show that quantum
effects indeed change collective universal behavior.

Introduction.— In reaction-diffusion (RD) models classical reactants, or
particles, are transported by diffusion and react when they meet, see, e.g.,
Refs. Vladimir (1997); Hinrichsen (2000); Henkel  _et al._ (2008). These are
paradigmatic non-equilibrium systems displaying universal dynamical properties
and stationary-state transitions from fluctuating phases to absorbing states,
i.e., states that once reached cannot be left. In one dimension, in
particular, spatial fluctuations of the particle number dominate the kinetics
and both exact analytical results Toussaint and Wilczek (1983); Spouge (1988);
Privman (1994); Torney and McConnell (1983); Rácz (1985); Takayasu  _et al._
(1988) and dynamical field-theory renormalization calculations Doi (1976);
Täuber  _et al._ (2005); Mattis and Glasser (1998); Peliti (1986, 1985);
Täuber (2014); Grassberger and Scheunert (1980) have shown that the dynamical
critical behavior is universal and it is not captured by the mean-field
approximation. This is especially true in the diffusion-limited regime, i.e.,
when the diffusive mixing of the particles is not too strong Ovchinnikov and
Zeldovich (1978); Kang and Redner (1984a); Toussaint and Wilczek (1983); Kang
and Redner (1984b); Privman (1994). In the opposite reaction-limited regime,
where the diffusive motion is fast, the density of reactants rapidly
uniformize (leading to the alternative name of well-stirred-mixture
approximation) and one recovers mean-field results Hinrichsen (2000); Vladimir
(1997); Täuber  _et al._ (2005); Kang and Redner (1985); Privman and Grynberg
(1992).

Quantum effects can alter the universal properties of absorbing-state phase
transitions. This has been shown for Markovian open quantum systems Griessner
_et al._ (2006); Diehl  _et al._ (2008); Kraus  _et al._ (2008); Diehl  _et
al._ (2011); Tomadin  _et al._ (2011); Bardyn  _et al._ (2013); Pérez-
Espigares  _et al._ (2017); Buča  _et al._ (2020); Carollo and Lesanovsky
(2022), for systems with kinetic constraints Lesanovsky and Garrahan (2013);
Olmos  _et al._ (2014); Everest  _et al._ (2016); Marcuzzi  _et al._ (2016);
Buchhold  _et al._ (2017); Gutiérrez  _et al._ (2017a); Roscher  _et al._
(2018); Carollo  _et al._ (2019); Gillman  _et al._ (2019, 2020); Wintermantel
_et al._ (2020); Helmrich  _et al._ (2020); Nigmatullin  _et al._ (2021);
Kazemi and Weimer (2021); Carollo  _et al._ (2022) and for the quantum contact
process Carollo  _et al._ (2019); Jo  _et al._ (2021). Quantum dissipative RD
spin chains, where the diffusive motion is replaced by coherent hopping, have
been investigated in Ref. van Horssen and Garrahan (2015). However, results in
this and other works are limited to small systems, due to the complexity of
the numerical simulation of many-body quantum dynamics. As a consequence, very
little is known about the impact of quantum effects on universal aspects of RD
dynamics and on absorbing-state phase transitions.

In this manuscript, we make progress in this direction, deriving exact
analytical results for the case of reaction-limited open quantum RD processes
in fermionic chains. We consider a series of prototypical reaction processes,
such as annihilation $A+A\to\emptyset$, coagulation $A+A\to A$, and branching
$A\to A+A$ (see Fig. 1), and show that the reaction-limited regime of quantum
RD models cannot be described within a mean-field approach, in stark contrast
to the classical settings. We demonstrate that the presence of quantum effects
strongly affects the approach to stationarity and the stationary state itself.
For annihilation and coagulation, the density of particles features an
algebraic (power-law) decay. This power law changes and may deviate from the
mean-field predictions when the initial state of the dynamics features quantum
coherence. In the presence of the branching process, quantum RD models display
an absorbing-state phase transition. Here, annihilation processes that couple
to coherent superpositions of adjacent particle pairs lead to the emergence of
dark states which are locally protected against dissipation. These local dark
states, which are not captured by the mean-field approach, establish quantum
correlations between fermionic particles.

Our analysis is performed by exploiting the time-dependent generalized Gibbs
ensemble method (TGGE) Lange  _et al._ (2018); Mallayya  _et al._ (2019);
Lange  _et al._ (2017); Lenarčič  _et al._ (2018), which naturally leads to
large-scale Boltzmann-like equations. The latter provides an exact description
for the reaction-limited regime in the thermodynamic limit. Our analytical
findings show that quantum effects lead to rich non-equilibrium behavior,
significantly different from that of classical systems. Our results connect to
the physics of cold atoms, where losses are of central experimental Syassen
_et al._ (2008); Burrows  _et al._ (2017); Bouchoule and Schemmer (2020);
Traverso  _et al._ (2009); Yamaguchi  _et al._ (2008); Kinoshita  _et al._
(2005); Söding  _et al._ (1999); Tolra  _et al._ (2004) and theoretical
García-Ripoll  _et al._ (2009); Everest  _et al._ (2014); Rossini  _et al._
(2021); Rosso  _et al._ (2022); Bouchoule  _et al._ (2020); Bouchoule and
Dubail (2021, 2022); Rosso  _et al._ (2023) relevance.

Quantum reaction-diffusion models— We consider fermionic quantum chains of
length $L$. Each site $j$ can be either occupied
$n_{j}\ket{\cdots\bullet_{j}\cdots}=\ket{\cdots\bullet_{j}\cdots}$ or empty
$n_{j}\ket{\cdots\circ_{j}\cdots}=0$, where $n_{j}=c_{j}^{\dagger}c_{j}$ and
the operators $c_{j},c_{j}^{\dagger}$ obey the fermionic anticommutation
relations $\\{c_{j},c_{j^{\prime}}^{\dagger}\\}=\delta_{j,j^{\prime}}$. The
fermionic statistics prevents double occupancy of lattice sites, typically
assumed in RD classical models Hinrichsen (2000); Henkel  _et al._ (2008);
Vladimir (1997). The dynamics is ruled by the quantum master equation Gorini
_et al._ (1976); Lindblad (1976); Breuer and Petruccione (2002) ($\hbar=1$
henceforth)

$\dot{\rho}(t)=-i[H,\rho(t)]+\mathcal{D}[\rho(t)].$ (1)

Here, we assume that the diffusive motion of the particles in classical RD
models is replaced by coherent hopping, which is accounted for by the quantum
Hamiltonian

$H=-\Omega\sum_{j=1}^{L}(c_{j}^{\dagger}c_{j+1}+c_{j+1}^{\dagger}c_{j})\,,$
(2)

with $\Omega$ the hopping rate [cf. Fig. 1(a)]. Such Hamiltonian is
diagonalized with Fourier-space fermionic operators
$\hat{c}_{k},\hat{c}_{k}^{\dagger}$, where $k$ is the quasi-momentum, and the
number operators $\hat{n}_{k}=\hat{c}_{k}^{\dagger}\hat{c}_{k}$ Franchini
(2017). It conserves the total number $N=\sum_{j}n_{j}=\sum_{k}\hat{n}_{k}$ of
particles: $[H,N]=0$. The irreversible reaction processes are encoded in the
dissipator $\mathcal{D}$. It takes the (Lindblad) form Gorini  _et al._
(1976); Lindblad (1976); Breuer and Petruccione (2002)

$\mathcal{D}[\rho]=\sum_{j,\nu}\left[L_{j}^{\nu}\rho{L_{j}^{\nu}}^{\dagger}-\frac{1}{2}\left\\{{L_{j}^{\nu}}^{\dagger}L_{j}^{\nu},\rho\right\\}\right],$
(3)

where $L_{j}^{\nu}$ are local jump operators. We consider four different
reactions, labelled by the parameter $\nu$, which are sketched in Fig. 1(a).
The first is binary annihilation, $A+A\to\emptyset$, of a pair of neighboring
particles (rate $\Gamma_{\alpha}$), which is described by the jump operators

$L_{j}^{\alpha}=L_{j}^{\alpha}(\theta)=\sqrt{\Gamma_{\alpha}}c_{j}(\cos\theta\,c_{j+1}-\sin\theta\,c_{j-1}).$
(4)

The sum of the two terms, whose balance is controlled by the angle
$\theta\in[0,\pi)$, allows for the possibility that interference between two
quantum mechanical amplitudes contributes to the pair annihilation process.
Such structure naturally emerges in the Bose-Hubbard model subject to strong
two-body losses. In this limit, the model can be mapped to free fermions (2)
with weak, $\Gamma_{\alpha}\ll\Omega$, two-body losses (4), as shown in Refs.
Syassen  _et al._ (2008); García-Ripoll  _et al._ (2009); Rossini  _et al._
(2021). The classical-incoherent annihilation process is recovered for
$\theta=0,\pi/2$. The second reaction is coagulation, $A+A\to A$, of a
particle upon meeting a neighbouring one (rate $\Gamma_{\gamma}/2$), with jump
operators

$L_{j}^{\gamma\pm}=\sqrt{\Gamma_{\gamma}/2}\,c_{j}n_{j\pm 1}.$ (5)

The third reaction is one-body annihilation, $A\to\emptyset$, (rate
$\Gamma_{\delta}$) with jump operators

$L_{j}^{\delta}=\sqrt{\Gamma_{\delta}}\,c_{j}\,.$ (6)

Figure 1: Quantum RD dynamics in the reaction limited regime. (a) Quantum
chain with sites that can either be occupied by a fermion,
$\ket{\cdots\bullet_{j}\cdots}$, or empty $\ket{\cdots\circ_{j}\cdots}$.
Particles can hop between nearest-neighboring sites with hopping rate
$\Omega$, Eq. (2). Dissipation consists of irreversible reactions at rate
$\Gamma_{\nu}$, Eqs. (4)-(6). The parameter $\theta$ controls coherent
superposition from pair annihilation events. (b) In the reaction-limited
regime, $\Gamma\ll\Omega$, reaction dynamics is slow and takes place on the
timescale $\thicksim\Gamma^{-1}$. Fast hopping rapidly smooths out spatial
fluctuations (highlighted in red), due to local reactions, and the state of
the systems is described by a homogeneous GGE($\tau$) (blue horizontal lines)
at any rescaled time $\tau=\Gamma t$. (c) The total particle density
$\braket{n}_{\mathrm{GGE}}(\tau)$ decays algebraically in rescaled time $\tau$
(blue points) for annihilation or coagulation with exponent dependent on
initial state coherence. When branching is included an absorbing-state phase
transition to an active, finite density of particles, state can occur. The
latter displays correlation when $\theta\neq 0,\pi/2$.

These three reactions break number conservation and, due to continued particle
loss, drive the system towards an absorbing state devoid of particles. To
establish a non-trivial steady state, we consider a fourth reaction, namely
branching, $A\to A+A$. This process allows for creation of a particle in the
neighborhood of an occupied site (rate $\Gamma_{\beta}/2$)

$L_{j}^{\beta\pm}=\sqrt{\Gamma_{\beta}/2}\,c_{j}^{\dagger}n_{j\pm 1}.$ (7)

The competition between the branching process and one-body annihilation (as in
the contact process Hinrichsen (2000); Henkel  _et al._ (2008)) gives rise to
a nonequilibrium absorbing-state phase transition, from the empty state to a
stationary active one with finite density of particles. Coagulation (5) and
branching (7) can be experimentally implemented in the facilitation regime
Lesanovsky and Garrahan (2014) of cold-atomic gases dressed with Rydberg
interactions Valado  _et al._ (2016); Gutiérrez  _et al._ (2017b);
Wintermantel  _et al._ (2021). For convenience, in the following when multiple
reactions are present, we rescale rates as $\Gamma_{\nu}=\Gamma\nu$, so that
$\Gamma$ sets the timescale of the dissipation, while the dimensionless
parameters $\alpha$, $\beta$, $\gamma$ and $\delta$ encode the relative
strength of the reactions [see Fig. 1(b)-(c)].

There are two important timescales in the dynamics: the reaction time
$\thicksim\Gamma^{-1}$, which gives the typical time needed for neighbouring
particles to react, and the hopping time (or diffusion time in classical RD
models) $\thicksim\Omega^{-1}$, which sets the timescale for two reacting
particles to meet. In classical settings Hinrichsen (2000); Vladimir (1997),
the dynamics qualitatively changes depending on the ratio $\Gamma/\Omega$. The
regime with $\Gamma/\Omega\gg 1$ is named diffusion limited as the propagation
of particles is the limiting factor for reactions to occur. In this regime,
spatial fluctuations are relevant and in one dimension the total particle
density $\braket{n}(t)=\braket{N}(t)/L$ decays algebraically as
$\braket{n}(t)\thicksim(\Omega t)^{-1/2}$ Toussaint and Wilczek (1983); Spouge
(1988); Privman (1994); Torney and McConnell (1983); Rácz (1985); Takayasu
_et al._ (1988); Ovchinnikov and Zeldovich (1978); Kang and Redner (1984a);
Toussaint and Wilczek (1983); Kang and Redner (1984b); Privman (1994), which
is slower than the corresponding mean field prediction
$\braket{n}_{\mathrm{MF}}(t)\thicksim(\Gamma t)^{-1}$ (note the different
rescaling of time).

The opposite regime, $\Gamma/\Omega\ll 1$, is the reaction-limited one. Here,
spatial fluctuations are irrelevant as fast motion makes the particle density
homogeneous in space. For classical systems Hinrichsen (2000); Vladimir
(1997); Täuber  _et al._ (2005); Kang and Redner (1985); Privman and Grynberg
(1992) this regime is described by law of mass action rate equations, which
assert that the rate of change of reactants is proportional to the product of
their global densities. This approach disregards spatial correlations among
particles and it indeed reproduces the mean-field result
$\braket{n}_{\mathrm{MF}}(t)\thicksim(\Gamma t)^{-1}$. In what follows, we
consider the quantum analogue of this regime, see Fig. 1(b)-(c). As we show,
this regime is much richer than its classical counterpart, as coherent effects
give rise to collective behavior and quantum correlations beyond mean field.

Reaction-limited TGGE— For our quantum RD models, the reaction limited regime
$\Gamma/\Omega\ll 1$ is equivalent to a weak dissipation limit, which can be
analyzed with the recently proposed time-dependent generalized Gibbs ensemble
(TGGE) of Refs. Lange  _et al._ (2018); Mallayya  _et al._ (2019); Lange  _et
al._ (2017); Lenarčič  _et al._ (2018). Due to fast hopping, one can consider
the state of the system $\rho(t)$ to be relaxed with respect to the stationary
manifold of the Hamiltonian, $[H,\rho(t)]=0$, at any time $t$. The dynamics of
$\rho(t)$ within this manifold is set by the timescale $\Gamma^{-1}$ and it is
determined by the dissipation. This aspect is pictorially shown in Fig. 1(b).
The TGGE approach then makes an ansatz among the set of relaxed states of the
Hamiltonian, which is the GGE, see, e.g., Refs. Essler and Fagotti (2016);
Vidmar and Rigol (2016). In the specific case of the Hamiltonian (2), the GGE
can be written as

$\rho_{\mathrm{GGE}}(t)=\frac{1}{\mathcal{Z}(t)}\mbox{exp}\left(-\sum_{k}\lambda_{k}(t)\hat{n}_{k}\right),$
(8)

where $\mathcal{Z}(t)=\prod_{k}[1+e^{-\lambda_{k}(t)}]$. The GGE state (8)
describes averages $\braket{\dots}_{\mathrm{GGE}}(t)$ of local observables in
the thermodynamic limit. It is entirely fixed from the knowledge of the
Lagrange multipliers $\lambda_{k}(t)$ or, equivalently, of the occupation
functions $\braket{\hat{n}_{q}}_{\mathrm{GGE}}(t)=C_{q}(t)$, which obey the
equations Rossini  _et al._ (2021); Rosso  _et al._ (2022); Bouchoule  _et
al._ (2020); Rosso  _et al._ (2023)

$\frac{\mbox{d}C_{q}(t)}{\mbox{d}t}=\sum_{j,\nu}\braket{{L_{j}^{\nu}}^{\dagger}[\hat{n}_{q},L_{j}^{\nu}]}_{\mathrm{GGE}}(t),\,\quad\forall
q.$ (9)

The solution $C_{q}(\tau)$ of this equation clearly depends on the rescaled
time $\tau=\Gamma t$, consistently with the above discussion on the reaction-
limited regime. The equation of motion (9) describes the large-scale dynamics
of the system and it has a structure akin to the Boltzmann equation. The right
hand side can be, crucially, exactly computed in the GGE state (8) through
Wick’s theorem. To explore the impact of quantum-coherent effects on the RD
dynamics, we consider two different initial conditions for Eq. (9). The first
is the coherent Fermi-sea (FS) state with density-filling $0<n_{0}\leq 1$:
$C_{q}(t=0)=1$ if $q\in[-\pi n_{0},\pi n_{0}]$, and zero otherwise. The second
is the incoherent state $\rho_{0}=\mbox{exp}(-\lambda N)/\mathcal{Z}_{0}$,
with a flat initial distribution in momentum space, $C_{q}(0)=n_{0}$.

Figure 2: Dynamics and active phase in quantum reaction-limited RD systems.
(a) Log-log plot of the particle density $\braket{n}_{\mathrm{GGE}}(\tau)$ as
a function of the rescaled time $\tau=\Gamma_{\alpha}t$ for the binary
annihilation reaction (4) with $\theta=0$. In the top-blue curve, the initial
state is the coherent Fermi sea (FS) state with filling $n_{0}=0.7$. The
density decays asymptotically as a power law
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$. In the inset, the black
dashed curve is a power-law fit $\braket{n}_{\mathrm{GGE}}(\tau)=a\tau^{-b}$
performed over the time window $\tau\in[10^{6},10^{7}]$, with the resulting
fitting parameter for the exponent being $b=0.50025\pm 5\cdot 10^{-5}$. In the
red-dashed curve, the initial state is the incoherent state $\rho_{0}$ with
the same mean density $n_{0}=0.7$. In this case, the density is exactly
described by the mean-field (MF) law of mass action and
$\braket{n}_{\mathrm{GGE}}(\tau)=\braket{n}_{\mathrm{MF}}(\tau)\thicksim\tau^{-1}$.
(b) Log-log plot of the density of particles $\braket{n}_{\mathrm{GGE}}(\tau)$
as a function of $\tau=\Gamma_{\gamma}t$ for the coagulation reaction (5). The
top-blue curve corresponds to the FS initial state at filling $n_{0}=0.3$,
while the red-dashed one to the incoherent state $\rho_{0}$ at the same
filling. For the FS state, the asymptotic exponent
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1}$ is the same as in MF. (c)
Log-log plot of the density as a function of $\tau=\Gamma t$ for the CP with
pair annihilation Eqs. (4)-(6) and $\Gamma_{\gamma}=0$, from the FS initial
state at $n_{0}=0.7$. For $\beta>\delta$ an active stationary state is
reached. The associated stationary momentum distribution function
$C_{q}^{\mathrm{stat}}$ is shown in the inset as a function of $q$. (d)
Stationary correlations $G^{\mathrm{stat}}(2,\theta)$ at distance $2$ (left,
blue axis) and dark state contribution
$G^{\mathrm{dark}}(\theta)=\sin(2\theta)/2$ (right, red axis) in the CP as a
function of $\theta$. Parameters are $\beta=\alpha=1$, $\delta=0.5$.

Annihilation and coagulation— In Fig. 2(a), we plot, from Eq. (9) SM , the
particle density as a function of time for the pair annihilation reaction only
($\Gamma_{\gamma}=\Gamma_{\beta}=\Gamma_{\delta}=0$), Eq. (4) with $\theta=0$,
so that interference effects are excluded. The density decays as
$\braket{n}_{\mathrm{{GGE}}}(\tau=\Gamma_{\alpha}t)\thicksim(\Gamma_{\alpha}t)^{-1/2}$
for the FS initial state for any filling $n_{0}\neq 1$. The $1/2$ decay
exponent does not necessarily require considering pure states. It also occurs
for initial mixed states with an inhomogeneous in $q$ initial occupation
function $C_{q}(0)$ SM . In contrast, for the initial state $\rho_{0}$ and any
$n_{0}$, the law of mass action is recovered and the density is exactly given
by mean field,
$\braket{n}_{\mathrm{MF}}(\tau)\thicksim(\Gamma_{\alpha}t)^{-1}$. This shows
the relevance of coherent effects in the critical dynamics of the model, since
the algebraic decay of the density in the reaction-limited regime is not
described by the mean-field approximation whenever the initial state is
quantum coherent. In the latter case, the decay of the particle density is
slower than in the classical counterpart of the model, where only incoherent
initial states are possible and the long-time behavior of the density is
independent on the initial density $n_{0}$ Kroon  _et al._ (1993); Henkel  _et
al._ (1995, 1997).

In Fig. 2(b), we plot the particle density as a function of time for the
coagulation reaction only
($\Gamma_{\alpha}=\Gamma_{\beta}=\Gamma_{\delta}=0$), Eq. (5) SM . We find
that
$\braket{n}_{\mathrm{GGE}}(\tau=\Gamma_{\gamma}t)\thicksim(\Gamma_{\gamma}t)^{-1}$
both for the incoherent state $\rho_{0}$ and for the FS state. For all initial
conditions, we see mean-field like decay 111For the FS initial state, initial
coherences approximately rescale time by an $n_{0}$ dependent factor without
altering the asymptotic mean-field decay. , which is different from the
situation for pair annihilation at $\theta=0$, Fig. 2(a). This difference
between annihilation and coagulation processes is in stark contrast with
classical RD models, where both processes belong to the same universality
class and decay in the same way independently of initial conditions Kang and
Redner (1984a); Spouge (1988); Privman (1994); Henkel  _et al._ (1995); Krebs
_et al._ (1995); Simon (1995); Henkel  _et al._ (1997); Ben-Avraham and Brunet
(2005).

For quantum RD, only when starting from the incoherent initial state
$\rho_{0}$ annihilation and coagulation behave in a similar way. In fact, the
densities $\braket{n}_{\mathrm{GGE}}^{\mathrm{ann}}(\tau,n_{0})$ and
$\braket{n}_{\mathrm{GGE}}^{\mathrm{coag}}(\tau,n_{0})$ obey

$\braket{n}_{\mathrm{GGE}}^{\mathrm{coag}}(\tau,n_{0})=2\braket{n}_{\mathrm{GGE}}^{\mathrm{ann}}(\tau,n_{0}/2),$
(10)

for $\Gamma_{\alpha}=\Gamma_{\gamma}$. Equation (10) is proved noting that the
dynamics from the incoherent state $\rho_{0}$ according to Eq. (9) remains at
all times fully incoherent and the quantum master equation (1) can then be
mapped onto a classical master equation SM . For the coherent FS initial
state, off-diagonal elements of the density matrix $\rho(t)$ are relevant, the
quantum master equation does not reduce to its classical counterpart, and Eq.
(10) does not apply. This shows that the quantum RD annihilation and
coagulation processes do not generically belong to the same universality class
and they can display different asymptotic behavior.

Contact process— We now consider the contact process (CP) with pair
annihilation, cf. Eqs. (4)-(6) with $\Gamma_{\nu}=\Gamma\nu$
($\nu=\alpha,\beta,\delta$) and $\Gamma_{\gamma}=0$ and Fig. 1(c). In Fig.
2(c), we plot the density as a function of the rescaled time $\tau=\Gamma t$.
We find a phase transition between an absorbing and an active state: the
stationary-state density $\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}$ becomes
non-zero when $\beta>\beta_{c}$, with $\beta_{c}=\delta$ independent of
$\alpha$ and $\theta$. This $\beta_{c}$ is the same as that of the mean-field
classical CP Hinrichsen (2000); Henkel  _et al._ (2008). Furthermore, we find
that the associated critical exponents for the stationary density
$\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}\propto(\beta-\beta_{c})^{1}$ and
for the decay of the density at the critical point $\beta_{c}$,
$\braket{n}_{\mathrm{GGE}}\thicksim(\Gamma t)^{-1}$, are those of the (mean-
field) directed percolation universality.

Interestingly, however, the stationary state is strongly affected by the
quantum coherence introduced by the annihilation reaction in Eq. (4), beyond
what can be predicted by a mean-field approach. The inset of Fig. 2(c) shows
that the different quasi-momenta $q$ are not evenly populated in the
stationary state. This applies when $\theta\neq 0,\pi/2$. The non-trivial
structure of $C_{q}^{\mathrm{stat}}$ implies that the stationary state has
spatial correlations. To quantify this, we compute the two-point fermionic
correlation function
$G^{\mathrm{stat}}(x-y,\theta)=\braket{c^{\dagger}_{x}c_{y}}_{\mathrm{GGE}}^{\mathrm{stat}}$,
which for the mean-field (product) state would be zero unless $x=y$. We find
that $G^{\mathrm{stat}}(l,\theta)$ is non zero at even distances
$l=2,4,6\dots$ with a dominant contribution at $l=2$. The value of
$G^{\mathrm{stat}}(2,\theta)$ as a function of $\theta$ is shown in Fig. 2(d)
and is approximately equal to $A(\theta)=\varepsilon\sin(2\theta)/2$.
Considering only these dominant next-to-nearest-neighbor correlations, we can
identify the (approximate) Lagrange multipliers $\lambda_{q}^{\rm stat}$ for
the stationary GGE $\rho_{\mathrm{GGE}}^{\rm stat}$ expanding to first order
in $A(\theta)$ (since $\varepsilon$ is small as shown in Fig. 2(d)). One
obtains $\lambda^{\rm stat}_{q}=\lambda_{\rm MF}+\lambda_{2}\cos(2q)$ and
therefore $\rho_{\mathrm{GGE}}^{\rm stat}\propto e^{-\lambda_{\rm
MF}N-\lambda_{2}Q_{2}/2}$, with
$Q_{2}=\sum_{j}(c_{j}^{\dagger}c_{j+2}+c_{j+2}^{\dagger}c_{j})$. The
contribution $\lambda_{\rm MF}=\log(1/\braket{n}^{\rm stat}_{\mathrm{GGE}}-1)$
represents the mean-field component of the state, while
$\lambda_{2}=-A(\theta)/[\braket{n}^{\rm
stat}_{\mathrm{GGE}}(1-\braket{n}^{\rm stat}_{\mathrm{GGE}})]$ accounts for
deviations from it. We show SM  that $\rho_{\mathrm{GGE}}^{\rm stat}$ can be
written in terms of an incoherent state plus a coherent correction, where
projectors onto the local dark states

$\ket{\psi}^{\mathrm{dark},\circ/\bullet}_{j}=\pm\cos\theta\ket{\bullet(\circ/\bullet)_{j}\circ}+\sin\theta\ket{\circ(\circ/\bullet)_{j}\bullet},$
(11)

emerge out of the uncorrelated mean-field state. The states
$\ket{\psi}^{\mathrm{dark},\circ/\bullet}_{j}$ are both dark with respect to
the annihilation process (4) centered in $j$, i.e.,
$L_{j}^{\alpha}(\theta)\ket{\psi}^{\mathrm{dark},\circ/\bullet}_{j}=0$.
Moreover, $\ket{\psi}^{\mathrm{dark},\bullet}_{j}$ is dark to branching (7) in
$j$ and is connected through one-body annihilation (6) in $j$ to the state
$\ket{\psi}^{\mathrm{dark},\circ}_{j}$. These local dark states determine the
correlations $G^{\mathrm{stat}}(2,\theta)$ in $\rho_{\mathrm{GGE}}^{\rm
stat}$, as shown in Fig. 2(d).

Summary— We provided a fully analytical treatment of quantum many-body RD
systems in their reaction-limited regime, where the irreversible reaction
rates are much smaller than the coherent hopping rate. While for classical RD
models this regime is well described by a mean-field approach, we have shown
that quantum RD displays instead much richer behaviour. In particular, for
annihilation, quantum coherence in the initial state can give rise to an
algebraic density decay whose power-law exponent differs from the mean-field
one. Furthermore, we have shown that quantum annihilation and coagulation do
not belong to the same universality class. For the contact process plus pair
annihilation, we have found that the stationary state can feature
correlations, which emerge as a consequence of destructive interference. This
inherently quantum feature gives rise to locally protected and correlated dark
states. The RD systems discussed here connect the soft-matter physics of
chemical reactions to that of cold atoms, where reactions translate into
dissipative particle losses or creations Burrows  _et al._ (2017); Bouchoule
and Schemmer (2020); Traverso  _et al._ (2009); Yamaguchi  _et al._ (2008);
Kinoshita  _et al._ (2005); Söding  _et al._ (1999); Tolra  _et al._ (2004);
García-Ripoll  _et al._ (2009); Everest  _et al._ (2014); Rossini  _et al._
(2021); Rosso  _et al._ (2022); Bouchoule  _et al._ (2020); Bouchoule and
Dubail (2021, 2022); Rosso  _et al._ (2023), which can be implemented via
Rydberg dressing Valado  _et al._ (2016); Gutiérrez  _et al._ (2017b);
Wintermantel  _et al._ (2021). Quantum reaction-diffusion systems are an ideal
benchmark to investigate the impact of quantum effects on large-scale
universal properties via numerical methods Carollo  _et al._ (2019); Gillman
_et al._ (2019); Jo  _et al._ (2021) and dynamical Keldysh-field theory
Kamenev (2023); Sieberer  _et al._ (2016).

Acknowledgements.— G.P. acknowledges support from the Alexander von Humboldt
Foundation through a Humboldt research fellowship for postdoctoral
researchers. We acknowledge financial support in part from EPSRC Grant No.
EP/R04421X/1, EPSRC Grant No. EP/V031201/1, and the Leverhulme Trust Grant No.
RPG-2018-181. We are also grateful for financing from the Baden-Württemberg
Stiftung through Project No. BWST_ISF2019-23 and for funding from the Deutsche
Forschungsgemeinsschaft (DFG, German Research Foundation) under Project No.
435696605, as well as through the Research Unit FOR 5413/1, Grant No.
465199066. F.C. is indebted to the Baden-Württemberg Stiftung for financial
support by the Eliteprogramme for Postdocs.

## References

  * Vladimir (1997) P. Vladimir,  _Nonequilibrium statistical mechanics in one dimension_ (Cambridge University Press, 1997).
  * Hinrichsen (2000) H. Hinrichsen, Adv. Phys. 49, 815 (2000).
  * Henkel  _et al._ (2008) M. Henkel, H. Hinrichsen, S. Lübeck, and M. Pleimling,  _Non-equilibrium phase transitions_, Vol. 1 (Springer, 2008).
  * Toussaint and Wilczek (1983) D. Toussaint and F. Wilczek, J. Chem. Phys. 78, 2642 (1983).
  * Spouge (1988) J. L. Spouge, Phys. Rev. Lett. 60, 871 (1988).
  * Privman (1994) V. Privman, Phys. Rev. E 50, 50 (1994).
  * Torney and McConnell (1983) D. C. Torney and H. M. McConnell, J. Phys. Chem. 87, 1941 (1983).
  * Rácz (1985) Z. Rácz, Phys. Rev. Lett. 55, 1707 (1985).
  * Takayasu  _et al._ (1988) H. Takayasu, I. Nishikawa, and H. Tasaki, Phys. Rev. A 37, 3110 (1988).
  * Doi (1976) M. Doi, J. Phys. A: Math. Gen. 9, 1479 (1976).
  * Täuber  _et al._ (2005) U. C. Täuber, M. Howard, and B. P. Vollmayr-Lee, J. Phys. A: Math. Gen. 38, R79 (2005).
  * Mattis and Glasser (1998) D. C. Mattis and M. L. Glasser, Rev. Mod. Phys. 70, 979 (1998).
  * Peliti (1986) L. Peliti, J. Phys. A: Math. Gen. 19, L365 (1986).
  * Peliti (1985) L. Peliti, J. Phys. France 46, 1469 (1985).
  * Täuber (2014) U. C. Täuber,  _Critical dynamics: a field theory approach to equilibrium and non-equilibrium scaling behavior_ (Cambridge University Press, 2014).
  * Grassberger and Scheunert (1980) P. Grassberger and M. Scheunert, Fortschritte der Phys. 28, 547 (1980).
  * Ovchinnikov and Zeldovich (1978) A. Ovchinnikov and Y. B. Zeldovich, J. Chem. Phys. 28, 215 (1978).
  * Kang and Redner (1984a) K. Kang and S. Redner, Phys. Rev. A 30, 2833 (1984a).
  * Kang and Redner (1984b) K. Kang and S. Redner, Phys. Rev. Lett. 52, 955 (1984b).
  * Kang and Redner (1985) K. Kang and S. Redner, Phys. Rev. A 32, 435 (1985).
  * Privman and Grynberg (1992) V. Privman and M. D. Grynberg, J. Phys. A. Math. Gen. 25, 6567 (1992).
  * Griessner  _et al._ (2006) A. Griessner, A. J. Daley, S. R. Clark, D. Jaksch, and P. Zoller, Phys. Rev. Lett. 97, 220403 (2006).
  * Diehl  _et al._ (2008) S. Diehl, A. Micheli, A. Kantian, B. Kraus, H. P. Büchler, and P. Zoller, Nat. Phys. 4, 878 (2008).
  * Kraus  _et al._ (2008) B. Kraus, H. P. Büchler, S. Diehl, A. Kantian, A. Micheli, and P. Zoller, Phys. Rev. A 78, 042307 (2008).
  * Diehl  _et al._ (2011) S. Diehl, E. Rico, M. A. Baranov, and P. Zoller, Nat. Phys. 7, 971 (2011).
  * Tomadin  _et al._ (2011) A. Tomadin, S. Diehl, and P. Zoller, Phys. Rev. A 83, 013611 (2011).
  * Bardyn  _et al._ (2013) C.-E. Bardyn, M. A. Baranov, C. V. Kraus, E. Rico, A. İmamoğlu, P. Zoller, and S. Diehl, New J. Phys. 15, 085001 (2013).
  * Pérez-Espigares  _et al._ (2017) C. Pérez-Espigares, M. Marcuzzi, R. Gutiérrez, and I. Lesanovsky, Phys. Rev. Lett. 119, 140401 (2017).
  * Buča  _et al._ (2020) B. Buča, C. Booker, M. Medenjak, and D. Jaksch, New J. Phys. 22, 123040 (2020).
  * Carollo and Lesanovsky (2022) F. Carollo and I. Lesanovsky, Phys. Rev. Lett. 128, 040603 (2022).
  * Lesanovsky and Garrahan (2013) I. Lesanovsky and J. P. Garrahan, Phys. Rev. Lett. 111, 215305 (2013).
  * Olmos  _et al._ (2014) B. Olmos, I. Lesanovsky, and J. P. Garrahan, Phys. Rev. E 90, 042147 (2014).
  * Everest  _et al._ (2016) B. Everest, M. Marcuzzi, J. P. Garrahan, and I. Lesanovsky, Phys. Rev. E 94, 052108 (2016).
  * Marcuzzi  _et al._ (2016) M. Marcuzzi, M. Buchhold, S. Diehl, and I. Lesanovsky, Phys. Rev. Lett. 116, 245701 (2016).
  * Buchhold  _et al._ (2017) M. Buchhold, B. Everest, M. Marcuzzi, I. Lesanovsky, and S. Diehl, Phys. Rev. B 95, 014308 (2017).
  * Gutiérrez  _et al._ (2017a) R. Gutiérrez, C. Simonelli, M. Archimi, F. Castellucci, E. Arimondo, D. Ciampini, M. Marcuzzi, I. Lesanovsky, and O. Morsch, Phys. Rev. A 96, 041602 (2017a).
  * Roscher  _et al._ (2018) D. Roscher, S. Diehl, and M. Buchhold, Phys. Rev. A 98, 062117 (2018).
  * Carollo  _et al._ (2019) F. Carollo, E. Gillman, H. Weimer, and I. Lesanovsky, Phys. Rev. Lett. 123, 100604 (2019).
  * Gillman  _et al._ (2019) E. Gillman, F. Carollo, and I. Lesanovsky, New J. Phys. 21, 093064 (2019).
  * Gillman  _et al._ (2020) E. Gillman, F. Carollo, and I. Lesanovsky, Phys. Rev. Lett. 125, 100403 (2020).
  * Wintermantel  _et al._ (2020) T. M. Wintermantel, Y. Wang, G. Lochead, S. Shevate, G. K. Brennen, and S. Whitlock, Phys. Rev. Lett. 124, 070503 (2020).
  * Helmrich  _et al._ (2020) S. Helmrich, A. Arias, G. Lochead, T. M. Wintermantel, M. Buchhold, S. Diehl, and S. Whitlock, Nature 577, 481 (2020).
  * Nigmatullin  _et al._ (2021) R. Nigmatullin, E. Wagner, and G. K. Brennen, Phys. Rev. Research 3, 043167 (2021).
  * Kazemi and Weimer (2021) J. Kazemi and H. Weimer, arXiv:211.05352 (2021).
  * Carollo  _et al._ (2022) F. Carollo, M. Gnann, G. Perfetto, and I. Lesanovsky, Phys. Rev. B 106, 094315 (2022).
  * Jo  _et al._ (2021) M. Jo, J. Lee, K. Choi, and B. Kahng, Phys. Rev. Research 3, 013238 (2021).
  * van Horssen and Garrahan (2015) M. van Horssen and J. P. Garrahan, Phys. Rev. E 91, 032132 (2015).
  * Lange  _et al._ (2018) F. Lange, Z. Lenarčič, and A. Rosch, Phys. Rev. B 97, 165138 (2018).
  * Mallayya  _et al._ (2019) K. Mallayya, M. Rigol, and W. De Roeck, Phys. Rev. X 9, 021027 (2019).
  * Lange  _et al._ (2017) F. Lange, Z. Lenarčič, and A. Rosch, Nat. Commun. 8, 1 (2017).
  * Lenarčič  _et al._ (2018) Z. Lenarčič, F. Lange, and A. Rosch, Phys. Rev. B 97, 024302 (2018).
  * Syassen  _et al._ (2008) N. Syassen, D. M. Bauer, M. Lettner, T. Volz, D. Dietze, J. J. García-Ripoll, J. I. Cirac, G. Rempe, and S. Dürr, Science 320, 1329 (2008).
  * Burrows  _et al._ (2017) K. A. Burrows, H. Perrin, and B. M. Garraway, Phys. Rev. A 96, 023429 (2017).
  * Bouchoule and Schemmer (2020) I. Bouchoule and M. Schemmer, SciPost Phys. 8, 60 (2020).
  * Traverso  _et al._ (2009) A. Traverso, R. Chakraborty, Y. N. Martinez de Escobar, P. G. Mickelson, S. B. Nagel, M. Yan, and T. C. Killian, Phys. Rev. A 79, 060702 (2009).
  * Yamaguchi  _et al._ (2008) A. Yamaguchi, S. Uetake, D. Hashimoto, J. M. Doyle, and Y. Takahashi, Phys. Rev. Lett. 101, 233002 (2008).
  * Kinoshita  _et al._ (2005) T. Kinoshita, T. Wenger, and D. S. Weiss, Phys. Rev. Lett. 95, 190406 (2005).
  * Söding  _et al._ (1999) J. Söding, D. Guéry-Odelin, P. Desbiolles, F. Chevy, H. Inamori, and J. Dalibard, Appl. Phys. B 69, 257 (1999).
  * Tolra  _et al._ (2004) B. L. Tolra, K. M. O’Hara, J. H. Huckans, W. D. Phillips, S. L. Rolston, and J. V. Porto, Phys. Rev. Lett. 92, 190401 (2004).
  * García-Ripoll  _et al._ (2009) J. J. García-Ripoll, S. Dürr, N. Syassen, D. M. Bauer, M. Lettner, G. Rempe, and J. I. Cirac, New J. Phys. 11, 013053 (2009).
  * Everest  _et al._ (2014) B. Everest, M. R. Hush, and I. Lesanovsky, Phys. Rev. B 90, 134306 (2014).
  * Rossini  _et al._ (2021) D. Rossini, A. Ghermaoui, M. B. Aguilera, R. Vatré, R. Bouganne, J. Beugnon, F. Gerbier, and L. Mazza, Phys. Rev. A 103, L060201 (2021).
  * Rosso  _et al._ (2022) L. Rosso, A. Biella, and L. Mazza, SciPost Phys. 12, 44 (2022).
  * Bouchoule  _et al._ (2020) I. Bouchoule, B. Doyon, and J. Dubail, SciPost Phys. 9, 44 (2020).
  * Bouchoule and Dubail (2021) I. Bouchoule and J. Dubail, Phys. Rev. Lett. 126, 160603 (2021).
  * Bouchoule and Dubail (2022) I. Bouchoule and J. Dubail, J. Stat. Mech.: Theory Exp. 2022, 014003 (2022).
  * Rosso  _et al._ (2023) L. Rosso, A. Biella, J. De Nardis, and L. Mazza, Phys. Rev. A 107, 013303 (2023).
  * Gorini  _et al._ (1976) V. Gorini, A. Kossakowski, and E. C. G. Sudarshan, J. Math. Phys. 17, 821 (1976).
  * Lindblad (1976) G. Lindblad, Comm. Math. Phys. 48, 119 (1976).
  * Breuer and Petruccione (2002) H.-P. Breuer and F. Petruccione,  _The theory of open quantum systems_ (Oxford University Press on Demand, 2002).
  * Franchini (2017) F. Franchini,  _An introduction to integrable techniques for one-dimensional quantum systems_ (Springer, 2017).
  * Lesanovsky and Garrahan (2014) I. Lesanovsky and J. P. Garrahan, Phys. Rev. A 90, 011603 (2014).
  * Valado  _et al._ (2016) M. M. Valado, C. Simonelli, M. D. Hoogerland, I. Lesanovsky, J. P. Garrahan, E. Arimondo, D. Ciampini, and O. Morsch, Phys. Rev. A 93, 040701 (2016).
  * Gutiérrez  _et al._ (2017b) R. Gutiérrez, C. Simonelli, M. Archimi, F. Castellucci, E. Arimondo, D. Ciampini, M. Marcuzzi, I. Lesanovsky, and O. Morsch, Phys. Rev. A 96, 041602 (2017b).
  * Wintermantel  _et al._ (2021) T. Wintermantel, M. Buchhold, S. Shevate, M. Morgado, Y. Wang, G. Lochead, S. Diehl, and S. Whitlock, Nat. Commun. 12, 103 (2021).
  * Essler and Fagotti (2016) F. H. Essler and M. Fagotti, J. Stat. Mech.: Theory Exp. 2016, 064002 (2016).
  * Vidmar and Rigol (2016) L. Vidmar and M. Rigol, J. Stat. Mech.: Theory Exp. 2016, 064007 (2016).
  * (78) See Supplemental Material at [URL will be inserted by publisher] for the details of the calculations.
  * Kroon  _et al._ (1993) R. Kroon, H. Fleurent, and R. Sprik, Phys. Rev. E 47, 2462 (1993).
  * Henkel  _et al._ (1995) M. Henkel, E. Orlandini, and G. Schutz, J. Phys. A: Math. Gen. 28, 6335 (1995).
  * Henkel  _et al._ (1997) M. Henkel, E. Orlandini, and J. Santos, Ann. Phys. 259, 163 (1997).
  * Note (1) For the FS initial state, initial coherences approximately rescale time by an $n_{0}$ dependent factor without altering the asymptotic mean-field decay.
  * Krebs  _et al._ (1995) K. Krebs, M. P. Pfannmüller, B. Wehefritz, and H. Hinrichsen, J. Stat. Phys. 78, 1429 (1995).
  * Simon (1995) H. Simon, J. Phys. A: Math. Gen. 28, 6585 (1995).
  * Ben-Avraham and Brunet (2005) D. Ben-Avraham and É. Brunet, J. Phys. A: Math. Gen. 38, 3247 (2005).
  * Kamenev (2023) A. Kamenev,  _Field theory of non-equilibrium systems_ (Cambridge University Press, 2023).
  * Sieberer  _et al._ (2016) L. M. Sieberer, M. Buchhold, and S. Diehl, Rep. Prog. Phys. 79, 096001 (2016).

SUPPLEMENTAL MATERIAL

Reaction-limited quantum reaction-diffusion dynamics in dissipative chains

Gabriele Perfetto1, Federico Carollo1, Juan P. Garrahan2,3, and Igor
Lesanovsky1,2,3

1Institut für Theoretische Physik, Universität Tübingen, Auf der Morgenstelle
14, 72076 Tübingen, Germany
2School of Physics, Astronomy, University of Nottingham, Nottingham, NG7 2RD,
UK.
3Centre for the Mathematics, Theoretical Physics of Quantum Non-Equilibrium
Systems, University of Nottingham, Nottingham, NG7 2RD, UK

In this Supplemental Material we provide details about the calculations
presented in the main text. In Sec. I, we discuss the TGGE ansatz for the free
fermionic hopping Hamilonian, given by Eqs. (8) and (2) of the main text,
respectively. In Sec. II, we specialize the discussion of the reaction-limited
TGGE dynamics to the annihilation, coagulation and contact process reactions.
In Sec. III, we eventually prove the mapping between the annihilation and the
coagulation dynamics in Eq. (10) of the main text.

##  I Reaction-limited TGGE ansatz for the fermion hopping Hamiltonian

We consider the free-fermionic hopping Hamiltonian in (2). We take henceforth
periodic boundary conditions $c_{j+L}=c_{j}$. This choice is without loss of
generality as we always consider the thermodynamic limit $L\to\infty$, where
the choice of boundary conditions does not matter. The Hamiltonian is
diagonalized by Fourier transform Franchini (2017)

$H=-2\Omega\sum_{k_{n}}\cos(k_{n})\hat{n}_{k_{n}}\qquad\mbox{with}\qquad\hat{n}_{k}=\hat{c}^{\dagger}_{k}\hat{c}_{k},$
(S1)

and the operators $\hat{c}_{k}$ in Fourier space defined as

$\hat{c}_{k_{n}}=\frac{1}{\sqrt{L}}\sum_{j=1}^{L}e^{-ik_{n}j}c_{j},\qquad\mbox{with
inverse}\qquad
c_{j}=\frac{1}{\sqrt{L}}\sum_{k_{n}}e^{ik_{n}j}\hat{c}_{k_{n}}.$ (S2)

Here $k_{n}=2\pi n/L$, with $n=1,2\dots L$, are the quasi-momenta. In the
remainder of this Supplemental Material, we denote summations over the quasi
momenta $\sum_{k_{n}}$ as $\sum_{k}$ for simplicity. The Hamiltonian in Eq.
(S1) clearly commutes with $\hat{n}_{k}$ for every $k$ value:
$[H,\hat{n}_{k}]=0$. The Hamiltonian is integrable and it possesses an
extensive number of conserved charges. The latter are linearly related to the
$\hat{n}_{k}$ operators, see, e.g., the discussion in Refs. Essler and
Fagotti, 2016; Vidmar and Rigol, 2016. The generalized-Gibbs ensemble
$\rho_{\mathrm{GGE}}$, describing the relaxation at long times under the
unitary dynamics of Eq. (S1), can be therefore written in terms of the
$\hat{n}_{k}$ as in Eq. (8) of the main text. In the reaction limited/weak
dissipation regime $\Gamma/\Omega\ll 1$, one promotes the GGE to be time
dependent $\rho_{\mathrm{GGE}}\to\rho_{\mathrm{GGE}}(t)$, as proposed in Refs.
Lange  _et al._ , 2018; Mallayya  _et al._ , 2019; Lange  _et al._ , 2017;
Lenarčič  _et al._ , 2018. It is then convenient to introduce the adimensional
time $\tau=\Gamma t$, in terms of which the TGGE ansatz is formulated as

$\lim_{\Gamma/\Omega\to
0}\rho(t=\tau/\Gamma)=\rho_{\mathrm{GGE}}(\tau)=\frac{1}{\mathcal{Z}(\tau)}\mbox{exp}\left(-\sum_{k}\lambda_{k}(\tau)\hat{n}_{k}\right),\qquad\mbox{with}\qquad\frac{\mathrm{d}\rho_{\mathrm{GGE}}(t)}{\mathrm{d}t}=\mathcal{D}[\rho_{\mathrm{GGE}}(t)].$
(S3)

The last equation follows from $[H,\rho_{\mathrm{GGE}}(t)]=0$. We emphasize
that the TGGE describes, in the thermodynamic limit $L\to\infty$, the slow
evolution taking place on the time scale $\Gamma^{-1}$ of the full quantum
state $\rho(t)$. Within this limit, the hopping time $\Omega^{-1}$ (diffusion
classically) is much smaller than the reaction time $\Gamma^{-1}$ so that the
reactants rapidly mix in space rendering an homogeneous locally in
(generalized) equilibrium state $\rho_{\mathrm{GGE}}(\tau)$. The state (S3) is
Gaussian and diagonal in momentum space. Its dynamics is therefore entirely
encoded in the two-point function
$C_{k}(\tau)=\braket{\hat{c}_{k}^{\dagger}\hat{c}_{q}}_{\mathrm{GGE}}(\tau)=\delta_{k,q}/(\mbox{exp}(\lambda_{q})+1)$.
In particular, from Eq. (S3), one has

$\frac{\mbox{d}C_{q}(t)}{\mbox{d}t}=\frac{1}{2}\sum_{j,\nu}\braket{{L_{j}^{\nu}}^{\dagger}[\hat{n}_{q},L_{j}^{\nu}]}_{\mathrm{GGE}}(t)+\braket{[{L_{j}^{\nu}}^{\dagger},\hat{n}_{q}]L_{j}^{\nu}}_{\mathrm{GGE}}(t)=\sum_{j,\nu}\braket{{L_{j}^{\nu}}^{\dagger}[\hat{n}_{q},L_{j}^{\nu}]}_{\mathrm{GGE}}(t),\quad\forall
q.$ (S4)

In the first equality we used the cyclic invariance of the trace, while in the
second equality that $[\hat{n}_{q},\rho_{\mathrm{GGE}}(t)]=0$. The above
equation is recognized as Eq. (9) of the main text. We notice that the
dissipation timescale $\Gamma$ appears as common factor on the right hand side
of Eq. (S4) and the solution $C_{q}(\tau)$ and the TGGE state
$\rho_{\mathrm{GGE}}(\tau)$ therefore depend on the rescaled time $\tau$, as
anticipated in the main text.

##  II Quantum reaction-limited dynamics

In this Section we specialize Eq. (9) of the main text to the various reaction
processes considered. In Subsec. II.1, we consider the binary annihilation
reaction (4). In Subsec. II.2, we consider the coagulation reaction (5). In
Subsec. II.3, we eventually consider the contact process with binary
annihilation in Eqs. (4), (6) and (7).

###  II.1 Annihilation

For the annihilation reaction $A+A\to\emptyset$ (4), we write the jump
operators $L_{j}^{\alpha}$ in Fourier space according to Eq. (S2) as

$L_{j}^{\alpha}(\theta)=\sqrt{\Gamma_{\alpha}}c_{j}(\cos\theta
c_{j+1}-\sin\theta
c_{j-1})=\sqrt{\Gamma_{\alpha}}\frac{1}{L}\sum_{k,k^{\prime}}e^{i(k+k^{\prime})j}(\cos\theta
e^{ik^{\prime}}-\sin\theta e^{-ik^{\prime}})\hat{c}_{k}\hat{c}_{k^{\prime}}.$
(S5)

The following commutation relation is then useful in the evaluation of the
commutator in Eq. (9)

$[\hat{n}_{q},\hat{c}_{k}\hat{c}_{k^{\prime}}]=-\hat{c}_{k}\hat{c}_{k^{\prime}}(\delta_{k,q}+\delta_{k^{\prime},q}).$
(S6)

Inserting Eqs. (S5) and (S6) into Eq. (S4) one gets

$\frac{\mbox{d}C_{q}(t)}{\mbox{d}t}=-\frac{\Gamma_{\alpha}}{L}\sum_{k,k^{\prime}}\braket{\hat{c}_{k^{\prime}}^{\dagger}\hat{c}_{k}^{\dagger}\hat{c}_{q}\hat{c}_{k+k^{\prime}-q}}_{\mathrm{GGE}}(t)\left(\cos\theta
e^{-ik^{\prime}}\\!-\\!\sin\theta
e^{ik^{\prime}}\right)\\!\\!\left(\cos\theta(e^{i(k+k^{\prime}-q)}-e^{iq})\\!-\sin\theta(e^{-i(k+k^{\prime}-q)}-e^{-iq})\right).$
(S7)

From the previous equation, it is clear that upon rescaling the time as
$\tau=\Gamma_{\alpha}t$, the solution $C_{q}(\tau)$ depends only on $\tau$.
The four-point fermionic function in the previous equation is evaluated
exploiting the fact that the TGGE state (S4) is Gaussian and therefore Wick
theorem applies:

$\braket{\hat{c}_{k^{\prime}}^{\dagger}\hat{c}_{k}^{\dagger}\hat{c}_{q}\hat{c}_{k+k^{\prime}-q}}_{\mathrm{GGE}}(\tau)=C_{k}(\tau)C_{k^{\prime}}(\tau)(\delta_{k,q}-\delta_{k^{\prime},q}),$
(S8)

leading to the equation

$\frac{\mbox{d}C_{q}(\tau)}{\mbox{d}\tau}=-\frac{1}{L}\sum_{k}g_{\theta}(k,q)C_{k}(\tau)C_{q}(\tau).$
(S9)

The function $g_{\theta}(k,q)$ is given by

$\displaystyle g_{\theta}(k,q)$
$\displaystyle=2(1-\cos(k-q))+\sin(2\theta)(2\cos(k+q)-\cos(2k)-\cos(2q)).$
(S10)

Equations (S9) and (S10) have been used with $\theta=0$ to produce the data in
Fig. 2(a) of the main text. We checked that the solution of Eqs. (S9) and
(S10) is stable upon increasing $L$ from $L=400,500$ and $600$ and therefore
that the thermodynamic limit is reached. Setting $\theta=0$ into the
expression for $g_{\theta}(k,q)$ one has for Eq. (S9) that

$\frac{\mbox{d}C_{q}(\tau)}{\mbox{d}\tau}=-2C_{q}(t)\braket{n}_{\mathrm{GGE}}(\tau)+\frac{2}{L}\sum_{k}\cos(k-q)C_{k}(\tau)C_{q}(\tau),$
(S11)

and for the density of reactants

$\braket{n}_{\mathrm{GGE}}(\tau)=\sum_{q}C_{q}(\tau)/L,$ (S12)

that

$\frac{\mbox{d}\braket{n}_{\mathrm{GGE}}(\tau)}{\mbox{d}\tau}=-2\braket{n}_{\mathrm{GGE}}^{2}(\tau)+\frac{2}{L^{2}}\sum_{k,k^{\prime}}\cos(k-k^{\prime})C_{k}(\tau)C_{k^{\prime}}(\tau).$
(S13)

The last equation is not closed for the density
$\braket{n}_{\mathrm{GGE}}(\tau)$ because of the presence of the second term
on the right hand side. The first term on the right hand side of (S13) is
exactly the mean-field law of mass action describing the reaction-limited
regime of classical annihilation RD dynamics Hinrichsen (2000); Täuber  _et
al._ (2005); Henkel  _et al._ (2008). The time integration of this
contribution simply yields

$\frac{\mbox{d}\braket{n}_{\mathrm{MF}}(\tau)}{\mbox{d}\tau}=-2\braket{n}_{\mathrm{MF}}^{2}(\tau)\to\braket{n}_{\mathrm{MF}}(\tau=\Gamma_{\alpha}t)=\frac{n_{0}}{1+2\Gamma_{\alpha}\,tn_{0}}.$
(S14)

The factor $2$ in the previous equation accounts for the fact that in each
annihilation reaction $2$ particles are lost. The function
$\braket{n}_{\mathrm{MF}}(\tau)$ is depicted in red-dashed in Fig. 2(a). The
second term on the right hand side of Eq. (S13) causes the departure shown in
Fig. 2(a) from the law of mass action prediction (S14). In particular, this
term is non-zero if and only if the momentum distribution function
$C_{q}(\tau)$ is not flat in $q$. This is achieved, for example, for the Fermi
sea initial state at density $n_{0}\neq 1$ and it causes the decay
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$ shown in the blue line
of Fig. 2(a). In the opposite case, where $C_{q}(\tau)$ is flat in the
momentum $q$, the second term in the right hand side of (S13) is identically
zero and the classical mean-field prediction $\braket{n}_{\mathrm{MF}}(\tau)$
(S14) is exactly retrieved. This is precisely what happens for the incoherent
initial state $\rho_{0}$, with $C_{q}(\tau=0)=n_{0}$ for any $q$.

We mention that the decay of the density in the quantum RD annihilation
dynamics ($\theta=0$) has been also studied in Ref. van Horssen and Garrahan,
2015 via numerical simulations of quantum-jump trajectories for system sizes
up to $L=22$. Therein, the fully occupied initial state is taken, $n_{0}=1$
with our notation, and the diffusion (hopping)-limited regime
$\Omega=\Gamma_{\alpha}$ is considered. The density is found to decay in this
limit algebraically as $\braket{n}(t)\thicksim t^{-b}$, with $1/2<b<1$. In
light of our results, we expect the exponent $b(\Omega/\Gamma_{\alpha})$ to
vary as a function of $\Omega/\Gamma_{\alpha}$ towards the value
$b_{\mathrm{MF}}=1$ in Eq. (S14) attained in the reaction-limited regime at
large $\Omega/\Gamma_{\alpha}$. Similar algebraic decays, with an exponent
varying with the Hamiltonian to dissipation strength $\Omega/\Gamma$, have
been numerically observed in Refs. Jo  _et al._ , 2021; Carollo  _et al._ ,
2022 for different types of kinetically-constrained open quantum dynamics.

In the case $\theta\neq 0,\pi/2$, i.e., away from the classical limit of the
annihilation reaction, one notices that quantum coherences are produced by the
reaction part of the dynamics itself. This translates into the fact that Eq.
(S9) produces a non-homogeneous momentum occupation function $C_{q}(\tau)$,
even in the case the initial distribution $C_{q}(0)$ is flat in $q$. As a
consequence of this, one observes a decay
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$ for every FS initial
state, even at unit filling $n_{0}=1$, and, more generally, also for the
initial state $\rho_{0}$ at arbitrary $n_{0}$. This observation is consistent
with the results derived in Refs. Rossini  _et al._ , 2021; Rosso  _et al._ ,
2022 where two-body atomic losses in one-dimensional bosonic gases in the
dissipative quantum Zeno regime have been addressed.

####  II.1.1 Annihilation dynamics from momentum-inhomogeneous GGE initial
states

We present here additional analyses and examples in order to further
corroborate the results of Subsec. II.1 concerning the annihilation decay
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$. Henceforth in this
subsection we take $\theta=0$ in Eq. (S5). We consider the case where the
initial state has the GGE form in Eq. (S3) with momentum inhomogeneous initial
Lagrange multipliers $\lambda_{k}(0)$ and momentum occupation function

$\rho_{\mathrm{GGE}}(0)=\frac{1}{\mathcal{Z}(0)}\mbox{exp}\left(-\sum_{k}\lambda_{k}(0)\hat{n}_{k}\right),\quad\mathcal{Z}(0)=\prod_{k_{n},n=1}^{L}\left(1+e^{-\lambda_{k}(0)}\right),\quad
C_{q}(0)=\frac{1}{\mathrm{exp}(\lambda_{q}(0))+1}.$ (S15)

We call states as in Eq. (S15) momentum-inhomogeneous GGE initial states as
they assume a GGE form and they allow for an initial occupation function
$C_{q}(0)$ not flat in $q$. This implies that the various quasi-momenta $q$
are not uniformly populated in the initial state. It is immediate to compute
the purity $\mathcal{P}=\mbox{Tr}[\rho^{2}]$, for states as in (S15), as

$\mathcal{P}_{\mathrm{GGE}}(0)=\prod_{k_{n},n=1}^{L}\mathcal{P}_{k_{n}}(0),\quad\mbox{with}\quad\mathcal{P}_{k}(0)=\frac{\mbox{Tr}(e^{-2\lambda_{k}(0)\hat{n}_{k}})}{\left(1+e^{-\lambda_{k}(0)}\right)^{2}}=\frac{1+e^{-2\lambda_{k}(0)}}{1+2e^{-\lambda_{k}(0)}+e^{-2\lambda_{k}(0)}}.$
(S16)

It is clear that $0\leq\mathcal{P}_{k}(0)\leq 1$, with the upper bound $1$
attained only if $\lambda_{k}(0)\to\pm\infty$. From the relation in Eq. (S15)
between $C_{q}(0)$ and $\lambda_{q}(0)$, this implies $C_{q}(0)\to 1$ when
$\lambda_{q}(0)\to-\infty$, and $C_{q}(0)\to 0$ when
$\lambda_{q}(0)\to+\infty$. Consequently, $\mathcal{P}_{\mathrm{GGE}}(0)<1$,
and the initial state is mixed, whenever the occupation function $C_{q}(0)$ is
such that $0<C_{q}(0)<1$ for at least one quasi-momentum $q$. This is the
case, for example, of the state $\rho_{0}=\mbox{exp}(\lambda
N)/\mathcal{Z}_{0}$, considered in the main text, having flat occupation
function: $C_{q}(0)=n_{0}<1$ for every $q$. Conversely, for pure states, with
$\mathcal{P}_{\mathrm{GGE}}(0)=1$, the occupation function can attain only the
values $0$ or $1$. The latter is precisely the limiting case of the Fermi-sea
initial pure state. Investigating states of the form (S15) therefore
generalizes the analysis of the main text by allowing to study mixed states
with a momentum inhomogeneous initial distribution, with the (pure) coherent
Fermi sea and the incoherent (momentum-homogeneous) $\rho_{0}$ states
retrieved as particular cases. In Fig. S1, we show the dynamics of the
momentum occupation function $C_{q}(\tau)$ and the particle density
$\braket{n}_{\mathrm{GGE}}(\tau)$ for the annihilation reaction, with
$\theta=0$, starting from two different mixed states (S15) with momentum
inhomogeneous initial distribution. In both cases, the density decays as
$\braket{n}_{\mathrm{GGE}}\thicksim\tau^{-1/2}$. This shows that the non mean-
field decay $\braket{n}_{\mathrm{GGE}}\thicksim\tau^{-1/2}$ does not
necessarily require purity equal one for the initial state (as it is the case
for the Fermi sea). Mixed initial states (S15) with purity
$\mathcal{P}_{\mathrm{GGE}}(0)<1$ give the decay
$\braket{n}_{\mathrm{GGE}}\thicksim\tau^{-1/2}$ as well, as long as the
initial momentum occupation function $C_{q}(0)$ is inhomogeneous in $q$. The
latter is the necessary condition for the decay exponent to be $1/2$, as
commented in Subsec. II.1 on the basis of Eqs. (S11)-(S14).

Figure S1: Annihilation dynamics for the occupation function and the particle
density. . (a) Momentum occupation function $C_{q}(\tau)$ as function of $q$
for increasing values of the rescaled time $\tau=\Gamma_{\alpha}t$ (from top
to bottom) for the annihilation reaction at $\theta=0$. The initial occupation
function is $C_{q}(0)=(1+\cos(q))/2$, with initial density $1/2$ (topmost
black curve). (b) Top-blue curve: log-log plot of the density
$\braket{n}_{\mathrm{GGE}}(\tau)$ associated to the occupation function
$C_{q}(\tau)$ in panel (a). The density decays asymptotically in time as
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$. The bottom red-dashed
line gives the mean-field prediction in Eq. (S14) with $n_{0}=1/2$ and it is
reported for comparison. (c) Momentum occupation function $C_{q}(\tau)$ as a
function of $q$ for increasing values of $\tau$ from the initial condition
$C_{q}(0)=0.7(1+\sin(q))/2$ (topmost black curve), and initial density
$n_{0}=0.35$. Note that in this case $C_{q}(0)$ is not invariant under quasi-
momentum reversal $q\to-q$ and therefore $C_{q}(\tau)$ evolves differently for
positive and negative values of $q$. (d) Top-blue curve: log-log plot of the
density $\braket{n}_{\mathrm{GGE}}(\tau)$ associated to $C_{q}(\tau)$ in (c).
The density decays asymptotically as
$\braket{n}_{\mathrm{GGE}}(\tau)\thicksim\tau^{-1/2}$ also in this case,
differently from the mean-field prediction (red-dashed curve).

###  II.2 Coagulation

We consider here the symmetric coagulation reaction (5), where the reactions
$A+A\to\emptyset+A$ (right coagulation) and $A+A\to A+\emptyset$ (left
coagulation) happen with the same rate $\Gamma_{\gamma}/2$. The generalization
of the following analysis to the asymmetric coagulation, where right and left
coagulation take place with different rates $\Gamma_{\gamma+}$ and
$\Gamma_{\gamma-}$, respectively, is straightforward and it does not alter
qualitatively our results.

The calculation for the symmetric coagulation process $L_{j}^{\gamma\pm}$ in
Eq. (5) is more involved than the one for annihilation as it involves three
fermion operators. The expression in Fourier space of the jump operators is

$L_{j}^{\gamma\pm}=\sqrt{\Gamma_{\gamma}/2}\,c_{j}n_{j\pm
1}=\sqrt{\Gamma_{\gamma}/2}\,c_{j}\,c_{j\pm 1}^{\dagger}c_{j\pm
1}=\sqrt{\frac{\Gamma_{\gamma}}{2}}\frac{1}{L^{3/2}}\sum_{k_{1},k_{2},k_{3}}e^{ik_{1}j}e^{-ik_{2}(j\pm
1)}e^{ik_{3}(j\pm 1)}\hat{c}_{k_{1}}\hat{c}_{k_{2}}^{\dagger}\hat{c}_{k_{3}}.$
(S17)

The commutator in Eq. (S4) can be again simplified using the following
commutation relations

$\displaystyle[\hat{n}_{q},\hat{c}_{k_{1}}]=-\delta_{k_{1},q}\hat{c}_{q}\hat{n}_{q}=-\delta_{k_{1},q}\hat{c}_{q},\quad\mbox{and}\quad[\hat{n}_{q},\hat{c}_{k_{2}}^{\dagger}\hat{c}_{k_{3}}]=\hat{c}_{k_{2}}^{\dagger}\hat{c}_{k_{3}}(\delta_{k_{2},q}-\delta_{k_{3},q}).$
(S18)

From Eqs. (S17) and (S18) one has

$\displaystyle[\hat{n}_{q},L_{j}^{\gamma\pm}]=\sqrt{\frac{\Gamma_{\gamma}}{2}}\frac{1}{L^{3/2}}$
$\displaystyle\left[\sum_{k,k^{\prime}}\hat{c}_{k}\hat{c}_{q}^{\dagger}\hat{c}_{k^{\prime}}e^{ikj}e^{-iq(j\pm
1)}e^{ik^{\prime}(j\pm 1)}-\sum_{k,k^{\prime}}e^{ikj}e^{-ik^{\prime}(j\pm
1)}e^{iq(j\pm 1)}\hat{c}_{k}\hat{c}_{k^{\prime}}^{\dagger}\hat{c}_{q}\right.$
$\displaystyle\left.-\sum_{k,k^{\prime}}e^{iqj}e^{-ik(j\pm
1)}e^{ik^{\prime}(j\pm
1)}\hat{c}_{q}\hat{c}^{\dagger}_{k}\hat{c}_{k^{\prime}}\right],$ (S19)

and therefore

$\displaystyle\sum_{j}(L_{j}^{\gamma\pm})^{\dagger}[\hat{n}_{q},L_{j}^{\gamma\pm}]$
$\displaystyle=\frac{\Gamma_{\gamma}}{2}\frac{1}{L^{3}}\sum_{j}\sum_{k_{1},k_{2},k_{3}}e^{-ik_{1}j}e^{ik_{2}(j\pm
1)}e^{-ik_{3}(j\pm
1)}\hat{c}^{\dagger}_{k_{3}}\hat{c}_{k_{2}}\hat{c}_{k_{1}}^{\dagger}\left[\sum_{k,k^{\prime}}\hat{c}_{k}\hat{c}_{q}^{\dagger}\hat{c}_{k^{\prime}}e^{ikj}e^{-iq(j\pm
1)}e^{ik^{\prime}(j\pm 1)}\right.$
$\displaystyle\left.-\sum_{k,k^{\prime}}e^{ikj}e^{-ik^{\prime}(j\pm
1)}e^{iq(j\pm
1)}\hat{c}_{k}\hat{c}_{k^{\prime}}^{\dagger}\hat{c}_{q}-\sum_{k,k^{\prime}}e^{iqj}e^{-ik(j\pm
1)}e^{ik^{\prime}(j\pm
1)}\hat{c}_{q}\hat{c}^{\dagger}_{k}\hat{c}_{k^{\prime}}\right]$
$\displaystyle=\frac{1}{L^{2}}\left[\sum_{k_{1},k_{3},k,k^{\prime}}e^{\pm
i(k_{1}-k)}(\hat{c}_{k_{3}}^{\dagger}\hat{c}_{k_{1}-k+k_{3}+q-k^{\prime}}\hat{c}_{k_{1}}^{\dagger}\hat{c}_{k}\hat{c}_{q}^{\dagger}\hat{c}_{k^{\prime}}-\hat{c}_{k_{3}}^{\dagger}\hat{c}_{k_{1}+k_{3}-k+k^{\prime}-q}\hat{c}_{k_{1}}^{\dagger}\hat{c}_{k}\hat{c}_{k^{\prime}}^{\dagger}\hat{c}_{q})\right]$
$\displaystyle-$
$\displaystyle\frac{1}{L^{2}}\left[\sum_{k_{1},k_{3},k,k^{\prime}}e^{\pm
i(k_{1}-q)}\hat{c}_{k_{3}}^{\dagger}\hat{c}_{k_{1}+k_{3}+k-k^{\prime}-q}\hat{c}_{k_{1}}^{\dagger}\hat{c}_{q}\hat{c}_{k}^{\dagger}\hat{c}_{k^{\prime}}\right].$
(S20)

One realizes from equation (S20), that in order to proceed further with the
calculation we need to compute six-point fermion correlation functions. This
is accomplished, similarly as in the case of the annhilation reaction dynamics
in Eq. (S8), exploiting the Gaussian structure of the TGGE state (S3) and
therefore the Wick theorem. We report the calculation for the first term in
the sum on the third line of Eq. (S20):

$\displaystyle\braket{\hat{c}_{k_{3}}^{\dagger}\hat{c}_{k_{1}-k+k_{3}+q-k^{\prime}}\hat{c}_{k_{1}}^{\dagger}\hat{c}_{k}\hat{c}_{q}^{\dagger}\hat{c}_{k^{\prime}}}_{\mathrm{GGE}}(\tau)$
$\displaystyle=C_{k_{1}}C_{q}C_{k_{3}}\delta_{k_{1},k}\delta_{q,k^{\prime}}+C_{k_{3}}C_{k_{1}}(1-C_{k})\delta_{k_{1},k^{\prime}}\delta_{k,q}+C_{k_{3}}C_{q}(1-C_{k_{1}})\delta_{k_{3},k}\delta_{k^{\prime},q}$
$\displaystyle-C_{k_{3}}C_{k_{1}}$
$\displaystyle(1-C_{q})\delta_{k_{3},k}\delta_{k_{1},k^{\prime}}+C_{k_{3}}(1-C_{q})(1-C_{k_{1}})\delta_{k_{3},k^{\prime}}\delta_{k,q}+C_{k_{3}}(1-C_{q})C_{k_{1}}\delta_{k_{3},k^{\prime}}\delta_{k_{1},k}.$
(S21)

The calculation of the TGGE expectation value of the other two terms on the
third line of Eq. (S20) works similarly and it is not reported for brevity. In
the previous equation, the argument of the momentum occupation function
$C_{q}(\tau)$ is the rescaled time $\tau=\Gamma_{\gamma}t$ and it is not
reported again for the sake of brevity. Taking the expectation value of Eq.
(S20) over the time-dependent GGE state, according to Eq. (S4), and using the
result (S21), after some algebra one obtains the equation

$\frac{\mbox{d}C_{q}(\tau)}{\mbox{d}\tau}=\braket{n}_{\mathrm{GGE}}^{2}(\tau)-2C_{q}(\tau)\braket{n}_{\mathrm{GGE}}(\tau)+\frac{2C_{q}(\tau)}{L}\sum_{k}\cos(q-k)C_{k}(\tau)-\frac{1}{L^{2}}\sum_{k,k^{\prime}}\cos(k-k^{\prime})C_{k}(\tau)C_{k^{\prime}}(\tau),$
(S22)

where $\tau=\Gamma_{\gamma}t$ and $\braket{n}_{\mathrm{GGE}}$ as in Eq. (S12).
The differential equation for the latter quantity can be written by summing
Eq. (S22) over all the modes $q$

$\frac{\mbox{d}\braket{n}_{\mathrm{GGE}}(\tau)}{\mbox{d}\tau}=-\braket{n}_{\mathrm{GGE}}^{2}(\tau)+\frac{1}{L^{2}}\left(\sum_{k,k^{\prime}}\cos(k-k^{\prime})C_{k}(\tau)C_{k^{\prime}}(\tau)\right).$
(S23)

Equation (S22) has been solved with $L=600$ to produce the data in Fig. 2(b)
of the main text. Equation (S23) has the very same structure as Eq. (S13) for
the annihilation reaction dynamics. The difference between the two reaction
processes lies, however, in the evolution equation for the momentum occupation
function $C_{q}(\tau)$, cf. Eq. (S9) with Eq. (S22). The classical reaction-
limited mean-field analysis is encoded into the first term on the right hand
side of (S23), which corresponds to the law of mass action for the coagulation
dynamics. The time integration of this contribution is

$\frac{\mbox{d}\braket{n}_{\mathrm{MF}}(\tau)}{\mbox{d}\tau}=-\braket{n}_{\mathrm{MF}}^{2}(\tau)\to\braket{n}_{\mathrm{MF}}(\tau=\Gamma_{\gamma}t)=\frac{n_{0}}{1+\Gamma_{\gamma}\,tn_{0}}.$
(S24)

We notice that in Eq. (S24) there is no factor $2$ (as in Eq. (S14) instead)
since each coagulation reaction depletes the number of particles by $1$. The
second term in Eq. (S23) is beyond the mean-field classical reaction-limited
description and it is not zero if and only if the momentum distribution is not
flat in momentum space. This is the case of the FS initial state, whose
dynamics is shown in blue in Fig. 2(b). In the opposite case, where
$C_{q}(\tau)$ is flat in momentum space, the mean-field solution (S24) is
retrieved. This is the case of the dynamics from the initial state $\rho_{0}$,
which we plot in red-dashed in Fig. 2(b).

###  II.3 Contact process with annihilation

In this Subsection we discuss the contact process with pair annihilation given
by Eqs. (4)-(6) without coagulation $\Gamma_{\gamma}=0$. We consider the case
of symmetric branching reactions (7), where $A+\emptyset\to A+A$ (right
branching) and $\emptyset+A\to A+A$ (left branching) happen with the same rate
$\Gamma_{\beta}/2$. The generalization to asymmetric branching is again
straightforward and it does not change qualitatively the results we are going
to present. We further rescale all the reactions’ rates
$\Gamma_{\nu}=\Gamma\nu$, with $\nu=\alpha,\beta$ and $\delta$, such that
$\Gamma$ sets the overall dissipation rate, while $\alpha$, $\beta$ and
$\delta$ encode the relative strength of the three reactions here considered.

The calculation for the branching jump operator $L_{j}^{\beta\pm}$ is similar
to the one explained in Subsec. II.2 for the coagulation dynamics. In
particular, one has in Fourier space

$(L_{j}^{\beta\pm})^{\dagger}=\sqrt{\frac{\Gamma\beta}{2}}\frac{1}{L^{3/2}}\sum_{k_{1},k_{2},k_{3}}e^{ik_{1}j}e^{-ik_{2}(j\pm
1)}e^{ik_{3}(j\pm 1)}\hat{c}_{k_{1}}\hat{c}^{\dagger}_{k_{2}}\hat{c}_{k_{3}}.$
(S25)

From Eqs. (S25) the calculation proceeds in a similar way as the one outlined
in Subsec. II.2. We report here just the final result for the sake of brevity

$\frac{\mbox{d}C_{q}(\tau)}{\mbox{d}\tau}=\beta\left(2\braket{n}_{\mathrm{GGE}}(\tau)-C_{q}(\tau)-\braket{n}_{\mathrm{GGE}}^{2}(\tau)+\frac{1}{L^{2}}\sum_{k,k^{\prime}}\cos(k-k^{\prime})C_{k}(\tau)C_{k^{\prime}}(\tau)\right)-\delta
C_{q}-\frac{\alpha}{L}\sum_{k}g_{\theta}(k,q)C_{k}(\tau)C_{q}(\tau),$ (S26)

with $\tau=\Gamma t$ and $\braket{n}_{\mathrm{GGE}}(\tau)$ in Eq. (S12).
Equation (S26) has been solved for $L=600$ to produce the data in Fig.
2(c)-(d). It is also instructive to look at the structure of the equation for
the density of particles:

$\frac{\mbox{d}\braket{n}_{\mathrm{GGE}}(\tau)}{\mbox{d}\tau}=\beta\braket{n}_{\mathrm{GGE}}(\tau)(1-\braket{n}_{\mathrm{GGE}}(\tau))-\delta\braket{n}_{\mathrm{GGE}}(\tau)+\frac{1}{L^{2}}\sum_{k,k^{\prime}}\left(\beta\cos(k-k^{\prime})-\alpha
g_{\theta}(k,k^{\prime})\right)C_{k}(\tau)C_{k^{\prime}}(\tau).$ (S27)

The part of the right hand side which solely depends on the density
$\braket{n}_{\mathrm{GGE}}(\tau)$ gives the result one would get within the
mean-field treatment of the classical reaction-limited RD dynamics. The terms
coupling different Fourier modes $C_{k}(\tau)$ and $C_{k^{\prime}}(\tau)$ go
beyond the latter description. From Eq. (S10) and (S27), the mean-field
prediction for the stationary density $n_{\mathrm{MF}}^{\mathrm{stat}}$ is
readily obtained

$\frac{\mbox{d}\braket{n}_{\mathrm{MF}}(\tau)}{\mbox{d}\tau}=-2\alpha\braket{n}_{\mathrm{MF}}^{2}(\tau)-\delta\braket{n}_{\mathrm{MF}}(\tau)+\beta\braket{n}_{\mathrm{MF}}(\tau)(1-\braket{n}_{\mathrm{MF}}(\tau))=0\to\braket{n}_{\mathrm{MF}}^{\mathrm{stat}}=\frac{\beta-\delta}{\beta+2\alpha},$
(S28)

which is defined ($\braket{n}_{\mathrm{MF}}^{\mathrm{stat}}\geq 0$) only if
$\beta>\beta_{c}=\delta$. In the reaction-limited regime, the critical point
of the absorbing-state phase transition therefore coincides with the one of
the classical CP (with branching (7) and decay (6) only) above its upper
critical dimension Hinrichsen (2000); Henkel  _et al._ (2008). The same
stationary state is furthermore obtained from the FS and the incoherent
initial state $\rho_{0}$ (the two initial states having the same density
$n_{0}$). The associated stationary density
$\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}$ is, however, strongly affected by
the coherences introduced by the annihilation reaction (4) at $\theta\neq
0(\pi/2)$. In particular, the stationary occupation function
$C_{q}^{\mathrm{stat}}=\lim_{\tau\to\infty}C_{q}(\tau)$, shown in the inset of
Fig. 2(c), is not flat as a function of $q$, which implies that the quantum
reaction-limited steady state displays spatial correlations beyond the
classical mean-field description. The stationary density achieved at long
times in the active phase, $\beta>\beta_{c}$, is consequently not given by Eq.
(S28) for $\theta\neq 0(\pi/2)$. In the case of Fig. 2(c), for example, we
find for $\theta=\pi/3$ (the other parameters are reported in the associated
caption) that $\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}\simeq 0.1706$, while
Eq. (S28) gives $\braket{n}^{\mathrm{stat}}_{\mathrm{MF}}=1/6\simeq 0.167$. In
order to investigate further the structure of the stationary state, we compute
the stationary correlation matrix
$G^{\mathrm{stat}}(x-y,\theta)=\lim_{\tau\to\infty}G(x-y,\theta,\tau)$, with:

$G(x-y,\theta,\tau)=\braket{c_{x}^{\dagger}c_{y}}_{\mathrm{GGE}}(\tau)=\frac{1}{L}\sum_{k,q}e^{iyq-
ixk}\braket{c_{k}^{\dagger}c_{q}}_{\mathrm{GGE}}(\tau)=\frac{1}{L}\sum_{q}e^{iq(y-x)}C_{q}(\tau).$
(S29)

The latter equation is nothing but the Fourier transform of $C_{q}(\tau)$.
Because of translational invariance, $G(x-y,\theta,\tau)$ is a function of the
distance $x-y$ between the two sites only. Moreover, since the initial
conditions investigated (the Fermi-sea and the incoherent state $\rho_{0}$)
and Eq. (S26) are invariant under quasi-momenta $q\to-q$ reversal,
$C_{q}(\tau)$ is at any time an even function of $q$. As a consequence, the
correlation matrix $G(l,\theta,\tau)$ is at any time $\tau$ real and symmetric
with respect to the origin:
$G^{\ast}(l,\theta,\tau)=G(l,\theta,\tau)=G(-l,\theta,\tau)$.

In Fig. S2, we plot $G^{\mathrm{stat}}(l,\theta)$ as a function of the
distance $l$ between the two sites, with the initial state taken as the FS at
filling $n_{0}=0.7$ (in the same way as in Fig. 2(c)-(d)). One can see that
$G^{\mathrm{stat}}(l,\theta)$ has a peak at $l=0$, whose magnitude corresponds
to the stationary density $\braket{n}^{\mathrm{stat}}_{\mathrm{GGE}}$ of the
active phase. In addition, $G^{\mathrm{stat}}(l,\theta)$ is non-zero at even
values of $l=\pm 2,\pm 4,\pm 6\dots$. This fact shows that the stationary
active state is not factorized in real space, as it would be within the mean-
field description of the classical reaction-limited dynamics. In the inset of
Fig. S2, we zoom in $G^{\mathrm{stat}}(l,\theta)$ away from $l=0$. The
dominant correlations clearly take place at distance $l=2$. Only in the case
$\theta=0\,(\pi/2)$, where the annihilation (4) reduces to its classical
limit, the steady state is uncorrelated and one recovers the classical mean-
field results:
$C_{q}^{\mathrm{stat}}=\braket{n}^{\mathrm{stat}}_{\mathrm{MF}}$ (flat in
momentum space) and $G^{\mathrm{stat}}(l,\theta)$ is zero for any $l\neq 0$.

Figure S2: Stationary correlations in the CP with binary annihilation. Plot of
$G^{\mathrm{stat}}(l,\theta)$ as a function of $l$ in the active stationary
state of the CP with pair annihilation (4) and without coagulation
$\Gamma_{\gamma}=0$ (cf. Eqs. (4)-(6)). The values of
$G^{\mathrm{stat}}(l,\theta)$ are represented by the red markers while the
blue dashed line is a guide for the eye. The value
$G^{\mathrm{stat}}(l=0,\theta)$ at $l=0$ corresponds to the stationary density
$\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}$, which is different from the mean-
field value $n_{\mathrm{MF}}^{\mathrm{stat}}$ as long as $\theta\neq 0,\pi/2$.
Fundamentally $G^{\mathrm{stat}}(l,\theta)$ is different from zero when
$\theta\neq 0,\pi/2$ at even distances $l=\pm 2,4,6,8,\dots$, showing that the
stationary state displays spatial correlations. In the inset, we zoom the
values of $G^{\mathrm{stat}}(l,\theta)$ in the interval $l\in[1,8]$, showing
that dominant correlations are at distance $l=\pm 2$, while correlations at
higher (even) distances are subleading. The parameters are analoguos to the
ones used in Fig. 2(c) of the main text: $\beta=\alpha=1$, $\delta=0.5$ and
$\theta=\pi/3$. The initial state is the FS at filling $n_{0}=0.7$. The same
result is obtained for the incoherent initial state $\rho_{0}$ at the same
density $n_{0}$.

We provide here the explicit expression of the GGE stationary state
$\rho_{\mathrm{GGE}}^{\mathrm{stat}}$ in order to better explain the relation
between the local dark states in Eq. (11) of the main text and the non-trivial
correlation function displayed in Fig. S2. As explained in the main text, the
stationary GGE can be written as $\rho_{\mathrm{GGE}}^{\mathrm{stat}}\propto
e^{-\lambda_{\mathrm{MF}}N-\lambda_{2}Q_{2}/2}$, with
$\lambda_{\mathrm{MF}}=\log(1/\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}-1)$
and
$\lambda_{2}=-\varepsilon\sin(\theta)/(2\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}(1-\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}))$.
We now expand the expression for $\rho_{\mathrm{GGE}}^{\mathrm{stat}}$ to
first order in $\varepsilon$ obtaining

$\rho_{\mathrm{GGE}}^{\mathrm{stat}}=\frac{e^{-\lambda_{\mathrm{MF}}N}}{\mathcal{Z}_{\mathrm{MF}}^{\mathrm{stat}}}-\frac{\lambda_{2}}{2}\frac{e^{-\lambda_{\mathrm{MF}}N}}{\mathcal{Z}_{\mathrm{MF}}^{\mathrm{stat}}}Q_{2}+\mathcal{O}(\varepsilon^{2}),\,\,\mbox{with}\,\,Q_{2}=\sum_{j}(c_{j}^{\dagger}c_{j+2}+c_{j+2}^{\dagger}c_{j}),\,\,\mbox{and}\,\,\mathcal{Z}^{\mathrm{stat}}_{\mathrm{MF}}=\mbox{Tr}(e^{-\lambda_{MF}N})=\prod_{j=1}^{L}(1+e^{-\lambda_{\mathrm{MF}}}).$
(S30)

In the previous step, we used the fact that $Q_{2}$ is a conserved charge of
the Hamiltonian and therefore $[Q_{2},N]=0$ and that
$e^{-\lambda_{\mathrm{MF}}N}$ is purely diagonal. The term
$e^{-\lambda_{\mathrm{MF}}N}Q_{2}$ is consequently purely off-diagonal so that
the normalization of $\rho_{\mathrm{GGE}}^{\mathrm{stat}}$ is
$Z_{\mathrm{MF}}^{\mathrm{stat}}$ to first order in $\varepsilon$. The first
term on the right hand side of the equation for
$\rho^{\mathrm{stat}}_{\mathrm{GGE}}$ is an incoherent mixture of states in
the fermionic Fock space spanned by $\ket{C}=\ket{C_{1},C_{2}\dots
C_{L}}=\ket{\circ_{1}\circ_{2}\dots\bullet_{L}}$, with
$N\ket{C}=N(C)\ket{C}=\sum_{j}N(C_{j})\ket{C}$, according to the factorized
probability measure $\propto e^{-\lambda_{MF}N(C)}$:

$\frac{e^{-\lambda_{\mathrm{MF}}N}}{\mathcal{Z}_{\mathrm{MF}}^{\mathrm{stat}}}=\sum_{C}\frac{e^{-\lambda_{\mathrm{MF}}N(C)}}{Z_{\mathrm{MF}}^{\mathrm{stat}}}\ket{C}\bra{C}=\sum_{C_{1},C_{2}\dots
C_{L}}\left(\prod_{j=1}^{L}P_{j}\right)\ket{C}\bra{C},\quad\mbox{and}\quad
P_{j}=\frac{e^{-\lambda_{\mathrm{MF}}N(C_{j})}}{1+e^{-\lambda_{\mathrm{MF}}}}.$
(S31)

One can eventually calculate the action of $Q_{2}$ onto the state (S31) which
leads to

$\displaystyle\rho^{\mathrm{stat}}_{\mathrm{GGE}}$
$\displaystyle=\frac{e^{-\lambda_{\mathrm{MF}}N}}{\mathcal{Z}_{\mathrm{MF}}^{\mathrm{stat}}}-\frac{\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}A(\theta)}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}(\ket{C_{<j-1},\circ_{j-1}\bullet_{j}\bullet_{j+1},C_{>j+1}}\bra{C_{<j-1},\bullet_{j-1}\bullet_{j}\circ_{j+1},C_{>j+1}}+\mbox{h.c.})$
$\displaystyle+\frac{(1-\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}})A(\theta)}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}(\ket{C_{<j-1},\circ_{j-1}\circ_{j}\bullet_{j+1},C_{>j+1}}\bra{C_{<j-1},\bullet_{j-1}\circ_{j}\circ_{j+1},C_{>j+1}}+\mbox{h.c.})\\!+\\!\mathcal{O}(\varepsilon^{2}),$
(S32)

with $A(\theta)=\varepsilon\sin(2\theta)/2$, as defined in the main text. In
the previous equation, we denoted with $C_{<j-1}$ ($C_{>j+1}$) the Fock state
preceding (following) the site $j-1$ ($j+1$), i.e., $\ket{C_{1},C_{2}\dots
C_{j-1}}$ ($\ket{C_{j+2},\dots C_{L}}$). We have also denoted with
$M_{j}=\prod_{l\neq j,j\pm 1}P_{l}$, the marginal distribution for all the
lattice sites but $j-1,j,j+1$. The relation between the previous equation and
the local dark states $\ket{\psi}_{j}^{\mathrm{dark},\circ/\bullet}$ in Eq.
(11) of the main text can be made more explicit upon rewriting Eq. (S32) as

$\displaystyle\rho^{\mathrm{stat}}_{\mathrm{GGE}}=\rho^{\mathrm{stat}}_{\mathrm{diag}}$
$\displaystyle+\frac{\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}\varepsilon}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\psi_{j}^{\mathrm{dark},\bullet},C_{>j+1}}\bra{C_{<j-1},\psi_{j}^{\mathrm{dark},\bullet},C_{>j+1}}$
$\displaystyle+\frac{(1-\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}})\varepsilon}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\psi_{j}^{\mathrm{dark},\circ},C_{>j+1}},\bra{C_{<j-1},\psi_{j}^{\mathrm{dark},\circ},C_{>j+1}}+\mathcal{O}(\varepsilon^{2}),$
(S33)

with

$\displaystyle\rho^{\mathrm{stat}}_{\mathrm{diag}}=$
$\displaystyle\frac{e^{-\lambda_{\mathrm{MF}}N}}{\mathcal{Z}_{\mathrm{MF}}^{\mathrm{stat}}}-\frac{\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}\varepsilon\sin^{2}\theta}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\circ_{j-1}\bullet_{j}\bullet_{j+1},C_{>j+1}}\bra{C_{<j-1},\circ_{j-1}\bullet_{j}\bullet_{j+1},C_{>j+1}}$
$\displaystyle-$
$\displaystyle\frac{\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}}\varepsilon\cos^{2}\theta}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\bullet_{j-1}\bullet_{j}\circ_{j+1},C_{>j+1}}\bra{C_{<j-1},\bullet_{j-1}\bullet_{j}\circ_{j+1},C_{>j+1}}$
$\displaystyle-$
$\displaystyle\frac{(1-\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}})\varepsilon\sin^{2}\theta}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\circ_{j-1}\circ_{j}\bullet_{j+1},C_{>j+1}}\bra{C_{<j-1},\circ_{j-1}\circ_{j}\bullet_{j+1},C_{>j+1}}$
$\displaystyle-$
$\displaystyle\frac{(1-\braket{n}_{\mathrm{GGE}}^{\mathrm{stat}})\varepsilon\cos^{2}\theta}{2}\sum_{j}\sum_{C_{<j-1},C_{>j+1}}\\!\\!\\!\\!\\!\\!\\!M_{j}\ket{C_{<j-1},\bullet_{j-1}\circ_{j}\circ_{j+1},C_{>j+1}}\bra{C_{<j-1},\bullet_{j-1}\circ_{j}\circ_{j+1},C_{>j+1}}.$
(S34)

The term $\rho_{\mathrm{diag}}^{\mathrm{stat}}$ is incoherent and gives zero
contribution to the correlation function,
$\mbox{Tr}[c_{x}^{\dagger}c_{y}\rho_{\mathrm{diag}}^{\mathrm{stat}}]=0$ for
$x\neq y$. The non-trivial correlations in Fig. 2(d) of the main text are
entirely determined by the second and third term in Eq. (S33) and, in
particular, by the coherences introduced by the projectors onto the dark
states $\ket{\psi}_{j}^{\mathrm{dark},\circ/\bullet}$ appearing therein. The
non-trivial structure of $\rho_{\mathrm{GGE}}^{\mathrm{stat}}$, determined by
the appearance of the conserved charge $Q_{2}$, is necessarily determined by
the dark states $\ket{\psi}_{j}^{\mathrm{dark},\circ/\bullet}$ of the
annihilation reaction. When $\theta=0,\pi/2$ and destructive interference in
Eq. (S5) is not possible, the dark states are not present and $Q_{2}$ is as
well absent in $\rho_{\mathrm{GGE}}^{\mathrm{stat}}$. The latter is in this
case solely determined by the conserved charge $N$ and it is trivially
factorized in space and uncorrelated.

##  III Mapping between annihilation and coagulation

In this Section we discuss for quantum reaction-limited RD systems the mapping
between annihilation (4) at $\theta=0$ (or, equivalently, $\pi/2$), and
coagulation (5). The mapping is valid for the incoherent initial state
$\rho_{0}$ and it is expressed by Eq. (10) of the main text, which relates the
density of reactants time evolution in the two reaction processes.

In this Section we use the Jordan-Wigner (JW) transformation to describe the
RD dynamics via spin operators Franchini (2017)

$c_{j}=S_{j}\sigma^{-}_{j},\quad
c_{j}^{\dagger}=S_{j}^{\dagger}\sigma^{+}_{j},\quad
S_{j}=\prod_{l=1}^{j-1}(-\sigma_{l}^{z}),\quad
n_{j}=c_{j}^{\dagger}c_{j}=\frac{1+\sigma^{z}_{j}}{2},\quad\sigma^{\pm}_{j}=\frac{\sigma_{j}^{x}\pm
i\sigma_{j}^{y}}{2},$ (S35)

with $\sigma^{x,y,z}_{j}$ the spin $1/2$ Pauli matrix at site $j$. One
realizes that the fermionic number operator
$n_{j}=c_{j}^{\dagger}c_{j}=\ket{\uparrow}_{j}\tensor{*}[_{j}]{\bra{\uparrow}}{}$
is identified with the projector onto the spin up state
$\sigma^{z}_{j}\ket{\uparrow}_{j}=+\ket{\uparrow}_{j}$. The Hermitian operator
$S_{j}=S_{j}^{\dagger}$ is usually named JW string. The annihilation reaction
(4) in terms of the spin operators reads as

$L_{j}^{\alpha}(\theta=0)=-\sqrt{\Gamma_{\alpha}}\sigma_{j}^{-}\sigma_{j+1}^{-},$
(S36)

while the coagulation reaction becomes (5)

$L_{j}^{\gamma\pm}=\sqrt{\Gamma_{\gamma}/2}\,S_{j}\sigma_{j}^{-}n_{j\pm 1}.$
(S37)

It is important to emphasize that Eq. (S37) contains the JW string $S_{j}$ and
it is therefore not local in the spin representation. We, however, show in
this Section that in the proof of Eq. (10) the string term $S_{j}$ in Eq.
(S37) does not matter. The Hamiltonian (2) with the JW transformation becomes
the $XX$ spin chain Franchini (2017)

$H=-\Omega\sum_{j=1}^{L}(\sigma_{j}^{-}\sigma_{j+1}^{+}+\sigma_{j}^{+}\sigma_{j+1}^{-}).$
(S38)

In order to prove Eq. (10), we introduce also jump operators $L_{j}^{D,R}$
($L_{j}^{D,L}$) giving incoherent hopping to the right (left)

$L_{j}^{D,R}=\sqrt{D}c_{j+1}^{\dagger}c_{j}=\sqrt{D}\sigma_{j+1}^{+}\sigma_{j}^{-},\quad
L_{j}^{D,L}=\sqrt{D}c_{j}^{\dagger}c_{j+1}=\sqrt{D}\sigma_{j}^{+}\sigma_{j+1}^{-},$
(S39)

at rate $D$. We remark that the boundary terms, $j=L$, in the Hamiltonian
(S38) and the boundary jump operators $L_{j=L}^{\alpha}(\theta=0)$,
$L_{j=L}^{\gamma+}$, $L_{j=1}^{\gamma-}$, $L_{j=L}^{D,R}$ and $L_{j=L}^{D,L}$
depend on the parity $(-1)^{N}$ of the fermionic number $N$. We do not write
these terms explicitly here, as the analysis of the reaction-limited regime
through the TGGE of Sec. I directly applies in the thermodynamic limit
$L\to\infty$. In this limit boundary terms can be neglected. In the proof of
Eq. (10), we consider the incoherent initial state $\rho_{0}$ with mean
density $n_{0}$:

$\rho_{0}=\frac{\mbox{exp}(-\lambda
N)}{\mathcal{Z}_{0}}=\prod_{j=1}^{L}\left(n_{0}n_{j}+(1-n_{0})(1-n_{j})\right)=\prod_{j=1}^{L}\left(n_{0}\ket{\uparrow}_{j}\tensor{*}[_{j}]{\bra{\uparrow}}{}+(1-n_{0})\ket{\downarrow}_{j}\tensor{*}[_{j}]{\bra{\downarrow}}{}\right).$
(S40)

The reaction-limited dynamics in Eq. (S3) from the initial state (S40) remains
incoherent at all times and diagonal in the classical basis spanned by product
states of the form, e.g.,
$\ket{C}=\ket{\uparrow\uparrow\downarrow\dots\uparrow}$. The reaction-limited
Lindblad dynamics (S3) can be therefore mapped to a classical master equation
by introducing the state vector $\ket{P(t)}$:

$\displaystyle\rho(t)$
$\displaystyle=\sum_{C}P_{C}(t)\ket{C}\bra{C}\to\ket{P(t)}=\sum_{C}P_{C}(t)\ket{C},$
(S41a) $\displaystyle\frac{\mathrm{d}\rho_{\mathrm{GGE}}(t)}{\mathrm{d}t}$
$\displaystyle=\mathcal{D}[\rho_{\mathrm{GGE}}(t)]\to\frac{\mathrm{d}P(C,t)}{\mathrm{d}t}=\sum_{C^{\prime}\neq
C}W(C^{\prime}\to
C)P(C^{\prime},t)-R(C)P(C,t)\to\frac{\mathrm{d}\ket{P(t)}}{\mathrm{d}t}=-\mathrm{H}\ket{P(t)}.$
(S41b)

Here, $\mathrm{H}$ is the Hamiltonian of the classical master equation (not to
be confused with the Hamiltonian $H$ ruling the original coherent dynamics
(1)-(3)) and it is given by

$\mathrm{H}=-\sum_{C}\sum_{C^{\prime}\neq C}W(C^{\prime}\to
C)\ket{C}\bra{C^{\prime}}+\sum_{C}R(C)\ket{C}\bra{C}.$ (S42)

Here, the transition $W(C^{\prime}\to C)$ and the escape rate $R(C)$ are
related to the jump operators in the dissipator $\mathcal{D}$ as

$W(C^{\prime}\to
C)=\sum_{j}|\braket{C}{L_{j}}{C^{\prime}}|^{2},\quad\mbox{and}\quad
R(C)=\sum_{C^{\prime}\neq C}W(C\to
C^{\prime})=\sum_{j}\braket{C}{L_{j}^{\dagger}L_{j}}{C}.$ (S43)

The initial state $\rho_{0}$ (S40) is mapped to the state $\ket{\rho_{0}}$

$\rho_{0}\to\ket{\rho_{0}}=\begin{pmatrix}n_{0}\\\
1-n_{0}\end{pmatrix}_{1}\otimes\begin{pmatrix}n_{0}\\\
1-n_{0}\end{pmatrix}_{2}\dots\otimes\begin{pmatrix}n_{0}\\\
1-n_{0}\end{pmatrix}_{L}.$ (S44)

We note that the string operator $S_{j}$ present in Eq. (S37) does not
contribute to the dynamics for a purely incoherent density matrix, as in Eq.
(S41), since $S_{j}^{2}=1$. For this reason, in the following, we do not
consider the JW string $S_{j}$ in the coagulation jump operators
$L_{j}^{\gamma\pm}$ (S37).

The mapping to the classical master equation (S41) applies both to the
annihilation (S36) $\mathcal{D}_{\alpha}$ and to the coagulation (S37)
$\mathcal{D}_{\gamma}$ dissipator. In both the cases, it can be shown that one
can include the incoherent hopping (S39) into the dissipators
$\mathcal{D}_{\alpha,D}$ $\mathcal{D}_{\gamma,D}$, as the reaction-limited
dynamics in Eqs. (S11) and (S22) for $C_{q}(\tau)$ for the incoherent
evolution (S41) is not changed upon including the jump operators (S39). The
advantage of doing this is that the dissipators $\mathcal{D}_{\alpha,D}$ and
$\mathcal{D}_{\gamma,D}$ map under Eq. (S41) to the classical Hamiltonians
$\mathrm{H}_{\alpha,D}^{\mathrm{ann}}$ and
$\mathrm{H}_{\gamma,D}^{\mathrm{coag}}$ of the corresponding classical
reaction diffusion systems Henkel  _et al._ (1995); Krebs  _et al._ (1995);
Simon (1995); Henkel  _et al._ (1997); Ben-Avraham and Brunet (2005). In the
latter case, the incoherent hopping accounts for the diffusive motion of the
classical reactants, with the rate $D$ in Eq. (S39) the diffusion constant.
For the classical annihilation-diffusion we have

$\mathrm{H}_{\alpha,D}^{\mathrm{ann}}=-\sum_{j}\left[D(\sigma_{j}^{-}\sigma_{j+1}^{+}+\sigma_{j}^{+}\sigma_{j+1}^{-})+\frac{\Delta_{\alpha}}{2}\sigma^{z}_{j}\sigma^{z}_{j+1}-\frac{\Gamma_{\alpha}}{4}(\sigma^{z}_{j}+\sigma^{z}_{j+1})-\frac{\Delta_{\alpha}-2D}{2}\right]-\Gamma_{\alpha}\sum_{j}\sigma_{j}^{-}\sigma_{j+1}^{-},$
(S45)

with $\Delta_{\alpha}=D-\Gamma_{\alpha}/2$. For the coagulation-diffusion
dynamics

$\mathrm{H}_{\gamma,D}^{\mathrm{coag}}=-\sum_{j}\left[D(\sigma_{j}^{-}\sigma_{j+1}^{+}+\sigma_{j}^{+}\sigma_{j+1}^{-})+\frac{\Delta_{\gamma}}{2}\sigma^{z}_{j}\sigma^{z}_{j+1}-\frac{\Gamma_{\gamma}}{4}(\sigma^{z}_{j}+\sigma^{z}_{j+1})-\frac{\Delta_{\gamma}-2D}{2}\right]-\Gamma_{\gamma}\sum_{j}(\sigma_{j}^{-}n_{j+1}+\sigma^{-}_{j+1}n_{j}),$
(S46)

and $\Delta_{\gamma}=D-\Gamma_{\gamma}/2$. At
$\Gamma_{\gamma}=\Gamma_{\alpha}$, the two Hamiltonians are related through a
similiarity transformation $B$, as shown in Ref. Krebs  _et al._ , 1995:

$\mathrm{H}_{\gamma=\alpha,D}^{\mathrm{coag}}=B\,\mathrm{H}_{\alpha,D}^{\mathrm{ann}}\,B^{-1},\quad\mbox{with}\quad
B=\bigotimes_{j=1}^{L}B_{j}\quad\mbox{and}\quad B_{j}=\begin{pmatrix}2&0\\\
-1&1\end{pmatrix}_{j},\quad B_{j}^{-1}=\begin{pmatrix}1/2&0\\\
1/2&1\end{pmatrix}_{j}.$ (S47)

The similarity matrix $B$ is built as the $L-$fold tensor product of the
matrix $B_{j}$ at the site $j$, which is the same for every lattice site. This
equation in the classical realm is considered as the hallmark of the
equivalence between annihilation and coagulation Hinrichsen (2000); Henkel
_et al._ (1995, 1997); Krebs  _et al._ (1995). Once Eq. (S47) is established,
the equivalence between the quantum reaction-limited annihilation dynamics and
the coagulation one is readily is established. Namely, one has

$\displaystyle\braket{n}_{\mathrm{GGE}}^{\mathrm{coag}}(\tau,n_{0})$
$\displaystyle=\mbox{Tr}[n_{j}\rho_{\mathrm{GGE}}^{\gamma}(t)]=\mbox{Tr}[n_{j}\rho_{\mathrm{GGE}}^{\gamma,D}(t)]=\braket{-}{n_{j}\,\mbox{exp}(-\mathrm{H}_{\gamma=\alpha,D}^{\mathrm{coag}}t)}{\rho_{0}}=\braket{-}{n_{j}\,B\,\mbox{exp}(-\mathrm{H}_{\alpha,D}^{\mathrm{ann}}t)\,B^{-1}}{\rho_{0}}$
$\displaystyle=2\braket{-}{n_{j}\,\mbox{exp}(-\mathrm{H}_{\alpha}^{\mathrm{ann}}t)}{\rho_{0}/2}=2\braket{n}_{\mathrm{GGE}}^{\mathrm{ann}}(\tau,n_{0}/2),$
(S48)

which corresponds to Eq. (10) of the main text and it implies that the density
of reactants decays with the same exponent in the two processes. Here,
$\ket{-}=\sum_{C}\ket{C}$ is the flat state, which implements the
normalization of the classical dynamics. In the second equality, we used the
fact that for an incoherent dynamics the introduction of the incoherent
hopping does not alter the quantum reaction-limited dynamics. In the third
equality, we used the mapping (S41) and, in the fourth equality, Eq. (S47). In
the fifth equality, we used that $n_{j}B_{j}=2n_{j}$, $\bra{-}B_{j}=\bra{-}$
and $B^{-1}\ket{\rho_{0}}=\ket{\rho_{0}/2}$. The derivation of Eq. (S48) can
be straightforwardly extended to density equal-time correlation functions
Henkel  _et al._ (1995); Krebs  _et al._ (1995); Henkel  _et al._ (1997); Ben-
Avraham and Brunet (2005).

Equations (S14) and (S24) do satisfy Eq. (S48). This derivation therefore
shows that the quantum reaction-limited dynamics of the annihilation (S36) and
coagulation (S37) processes from the initial state (S40) is exactly coincident
with the mean-field classical reaction-limited evolution. Departures from the
latter equivalence are of intrinsic quantum nature and they are determined
either by coherences in the initial state, as in the case of the FS initial
state, or by coherences introduced by the reactions, as in the case of Eq. (4)
at $\theta\neq 0,\pi/2$. In both these cases, the density matrix $\rho(t)$
develops in time coherences $\propto\ket{C}\bra{C^{\prime}}$ and the quantum
master equation (1) cannot be mapped to its classical counterpart, as in Eq.
(S41). In light of our results in Fig. 2(a)-(b), where the density of
particles decays with different asymptotic exponents in the two processes, we
conclude that Eq. (S48) does not hold when coherences are present and quantum
annihilation and coagulation generically display different asymptotic decays
(with different power-law exponents). This is in sharp contrast with the
classical case where annihilation and coagulation display analogous algebraic
decays irrespectively of the initial condition.

It is important to emphasize that the equivalence between annihilation and
coagulation is here proved in the reaction-limited regime in a weak way, i.e.,
in terms of the relation (S48) between the densities in the dynamics of the
two processes. In the classical case Hinrichsen (2000); Henkel  _et al._
(1995, 1997); Krebs  _et al._ (1995), the annihilation-coagulation equivalence
is proved in a stronger way in terms of the similarity relation (S47) between
the associated dynamical generators $H^{\mathrm{ann}}_{\alpha,D}$ and
$H^{\mathrm{coag}}_{\gamma=\alpha,D}$. Our results, however, do not rule out
the possibility of the equivalence between quantum annihilation and
coagulation, in the stronger sense that the associated Lindbladians (1)
(including therefore both the Hamiltonian $H$ in Eq. (2) and the dissipator
$\mathcal{D}$ (3)) are related through a similarity transformation. As a
matter of fact, in Ref. van Horssen and Garrahan, 2015, it has been suggested
that a similarity transformation relating the Lindbladians of the two process
does exist. This conclusion is drawn in Ref. van Horssen and Garrahan, 2015 on
the basis of exact-numerical diagonalization of the Lindbladian generators of
the two processes for up to $L=8$ spins. The existence of such similarity
relation between the two Lindblad generators does not, however, generically
imply a simple relation, as Eq. (S48), for the densities
$\braket{n}^{\mathrm{ann}}(t)$ and $\braket{n}^{\mathrm{coag}}(t)$ in the two
processes and therefore the same algebraic decay. One must, indeed, also
consider how the observable $n_{j}$ and the initial state
$B^{-1}\ket{\rho_{0}}$ transform under the similarity transformation $B$. In
order to establish Eq. (S48) it is, indeed, crucial that
$B^{-1}\ket{\rho_{0}}=\ket{\rho_{0}/2}$, i.e., the initial state is simply
transformed to a state of the same form but with halved density. Our results
seem to indicate that for the FS initial state the transformation rule is more
intricate, though a more in-depth analysis, which we leave for future studies,
is needed.