# Large $N$ Matrix Quantum Mechanics as a Quantum Memory

ChunJun Cao  Joint center for Quantum information and Computer science,
University of Maryland, College Park, MD 20740, USA  Institute for Quantum
Information and Matter, California Institute of Technology, Pasadena, CA
91125, USA  Gong Cheng  Maryland Center for Fundamental Physics, University of
Maryland, College Park, MD 20740, USA  Brian Swingle  Brandeis University,
Waltham, MA 02453, USA

###### Abstract

In this paper, we explore the possibility of building a quantum memory that is
robust to thermal noise using large $N$ matrix quantum mechanics models.
First, we investigate the gauged $SU(N)$ matrix harmonic oscillator and
different ways to encode quantum information in it. By calculating the mutual
information between the system and a reference which purifies the encoded
information, we identify a transition temperature, $T_{c}$, below which the
encoded quantum information is protected from thermal noise for a memory time
scaling as $N^{2}$. Conversely, for temperatures higher than $T_{c}$, the
information is quickly destroyed by thermal noise. Second, we relax the
requirement of gauge invariance and study a matrix harmonic oscillator model
with only global symmetry. Finally, we further relax even the symmetry
requirement and propose a model that consists of a large number $N^{2}$ of
qubits, with interactions derived from an approximate $SU(N)$ symmetry. In
both ungauged models, we find that the effects of gauging can be mimicked
using an energy penalty to give a similar result for the memory time. The
final qubit model also has the potential to be realized in the laboratory.

###### Contents

  1. 1 Introduction
    1. 1.1 Background on quantum codes
    2. 1.2 Background on quantum gauge theories
    3. 1.3 Organization and overview of results
  2. 2 Matrix model and set up
    1. 2.1 Matrix Quantum Mechanics
    2. 2.2 Physical states and the Knill-Laflamme Condition
    3. 2.3 Bath model
    4. 2.4 Mutual information diagnostic
  3. 3 Fully gauged model
    1. 3.1 Non-local model
    2. 3.2 Local model
  4. 4 Relaxing the gauge constraint
  5. 5 Large N spin model
  6. 6 Discussion and Conclusion
  7. A Proof of Theorem 1
    1. A.1 Early time approximation
    2. A.2 Finite time
  8. B aKL condition for single trace errors
    1. B.1 non-local model
    2. B.2 Local model
    3. B.3 Global symmetry
    4. B.4 Spin model
  9. C  Master equation in Lindblad form
  10. D Expectation value of the Casimir

##  1 Introduction

Quantum information is typically fragile when encountered in familiar physical
systems, thus giving rise to the intuition that quantum information is
inherently hard to preserve at macroscopic scales, e.g. in large systems or at
long times. The discovery of quantum error correction showed that this
intuition is incomplete: in the right kind of physical system, quantum
information can be retained for an arbitrarily long time. Maintaining the
information is typically an active process where we keep a vigil for errors
and periodically correct them. An error correcting code operating in this mode
forms an active quantum memory. Some codes may have the additional capability
to preserve the encoded information for a long time naturally, without any
active intervention. A code operating in this mode is called a passive or
self-correcting quantum memory—it is like a quantum hard drive. In this paper
and a forthcoming companion paper, we search for passive quantum memories in
the context of quantum gauge theories and models inspired by them.

Why gauge theories? For one thing, some of the most important currently known
codes can be viewed as gauge theories. For example, the toric code can be
viewed as a $Z_{2}$ lattice gauge theory with dynamical electric charges in
the extreme deconfined limit. More generally, deconfined lattice discrete
gauge theories are associated with a high degree of quantum entanglement and
topological order, both of which are linked with quantum error correction [1,
2, 3]. One feature of topological order is a robustness of the ground space to
perturbations of the Hamiltonian, which is suggestive of the potential for
passive error correction. And while the two-dimensional toric code does not
form a passive memory [4], a four-dimensional analog does [5]. Moreover,
phases of matter exhibiting topological order associated with a discrete non-
Abelian group can be used for fault tolerant quantum computation, a fact that
also hints at a degree of intrinsic protection [6].

From another perspective, quantum error correction has been linked with
quantum gravity in the context of holographic duality, also known as the
AdS/CFT correspondence [7]. The best understood examples of the correspondence
involve continuous non-Abelian gauge theories with large gauge groups, thus
providing another connection between non-Abelian gauge symmetry and quantum
error correction. In this context, one also finds that the boundary conformal
field theory with $SU(N)$ gauge symmetry can possess a low temperature
confined phase which is qualitatively distinct from the generic high
temperature physics in the deconfined phase. This boundary confinement-
deconfinement phase transition is dual to the Hawking-Page transition [8]
between a thermal AdS geometry and a AdS black hole geometry in the bulk [9].
It was suggested in [10] that when treating the low energy sector of AdS/CFT
as an approximate quantum error-correcting code under the influence of thermal
errors, the encoded information living in the bulk thermal AdS geometry can be
disturbed by the scattering of thermal particles with probability
$O(1/N^{2})$. However, in the high temperature phase, the encoded information
is scrambled by the black hole, thus inducing a large logical error. In the
large $N$ limit, the distinct phases correspond to the existence of a
threshold, such that logical information remains robust in the low temperature
phase where the errors occur less frequently111Note that unlike the tensor
network constructions [11] where one changes the system size by increasing the
number of layers, which roughly corresponds to changing the IR cutoff in the
bulk or the UV cutoff on the boundary theory when it is put on a lattice, the
size of the physical system in [10] is altered by $N$ of the $SU(N)$ symmetry
while keeping the number of “boundary sites” fixed. Heuristically, the
physical Hilbert space can be taken as placing $O(N^{2})$ qudits on each
boundary site..Therefore, it is natural to explore whether this type of
robustness also extends to other theories with such gauge symmetries.

With these motivations, we have undertaken an investigation of quantum error
correction in non-Abelian gauge theories with large gauge groups. We focus
here on $SU(N)$ groups with $N$ taken large. This large $N$ limit is like a
thermodynamic limit, and we expect that passive error correction can only
emerge in such a thermodynamic limit. We define and study quantum codes in a
variety of models, including models inspired by gauge theories but where we
relax the gauge constraint. We are seeking to answer two key questions: (1)
are there potentially useful lab-constructible quantum memories built from
gauge theories or derived models and (2) can naturally occurring systems, e.g.
QCD or a hypothetical holographic dual to quantum gravity in our universe,
naturally preserve quantum information for a long time? This paper focuses on
quantum mechanical systems with large gauge groups, while the companion paper
[12] will consider quantum field theories with large gauge groups.

The key results of the present paper are as follows. First, we discuss broadly
how a property called large $N$ factorization is related to an approximate
version of the Knill-Laflamme error correction condition [13, 14]. Assuming a
code space with the right properties can be defined and coupled to a thermal
bath satisfying reasonable assumptions, we then prove a theorem establishing
the existence of a memory lifetime which scales polynomially with $N$. We then
consider a sequence of models of increasing generality and show how the
conditions of the main theorem can be satisfied. For a certain fully gauged
model, we find that the conditions of the theorem are immediately satisfied;
for models where we do not gauge the symmetry or even make the symmetry
approximate, an energy penalty scaling like $\log N$ is needed to meet the
conditions of the main theorem. While here we focus on the construction of
such quantum memories and their memory times, there are other aspects of the
model that are relevant for practical implementation, such as whether quantum
information information can be efficiently encoded, recovered, and processed.
We will leave these discussions for future work.

In the remainder of the introduction, we first review the basic notions for
readers without the necessary background in quantum codes or quantum gauge
theories. We then give a high level overview of our results. We also note that
a few prior works have considered the role of large $N$ factorization and/or
gauge symmetry in defining quantum error correcting codes [15, 16]. There is
also a large literature on passive quantum memories from a variety of
approaches [17, 5, 2, 18]. See also [19] for a more comprehensive review.

###  1.1 Background on quantum codes

In this subsection, we review the key concepts of quantum error correction
needed for our discussion. This section is aimed at readers with a more high
energy background.

Mathematically, a quantum error correcting code is a subspace, the logical
space $\mathcal{C}$, of a larger physical Hilbert space $\mathcal{H}$. Errors
are modeled by a quantum process $\mathcal{E}$ and error correction is given
by another process $\mathcal{R}$ such that

$\mathcal{R}(\mathcal{E}(\rho))=\rho\quad\forall\rho\in\mathcal{C}.$ (1)

If the error channel is viewed as probabilistically applying errors from a set
$\\{E_{a}\\}$, then one way to formulate the necessary condition for the
existence of an error correction process is the Knill-Laflamme (KL) condition,

$\langle i|E_{a}^{\dagger}E_{b}|j\rangle=c_{ab}\delta_{ij},$ (2)

where $|i\rangle$ and $|j\rangle$ are any two elements of a basis for
$\cal{C}$. The intuition behind KL is that errors are unable to cause
transitions between distinct states in the code. We can also relax these exact
equalities to obtain instead an approximate quantum code.

In practice, such a code is instantiated within a concrete physical substrate,
with $\mathcal{H}$ being some model of the physical space of states of the
substrate. For example, a code ($\mathcal{C}$) might be implemented as a
subspace of the full energy levels ($\mathcal{H}$) of a superconducting
circuit (substrate). In some cases, the code space might correspond to the
ground space of some Hamiltonian for the substrate, but this is not required
in general.

A code is used as follows. Information is first encoded using some encoding
operation. This operation might correspond to a quantum circuit which takes a
fixed initial state in $\mathcal{H}$ and produces the desired code state in
$\mathcal{C}$. After encoding, the code is subject to errors, corresponding to
one application of the error channel $\mathcal{E}$. Then the recovery
procedure $\mathcal{R}$ is applied to recover the initially encoded
information. The condition above guarantees success provided the encoding
prepares a state in $\mathcal{C}$. Finally, the information can be read out by
a decoding operation. Note the distinction between recovery, which repairs the
information but keeps it encoded, and decoding, which reverses the encoding
map when applied to an error free state222Note that somewhat confusingly,
decoding in the context of active error correction can also refer to applying
some suitably chosen operations to the state based on syndrome measurements.
This process typically leave the information encoded, which corresponds to our
recovery operation in our current definition. .

To make these ideas more concrete, we can consider the simple goal of
information storage for a period of time (as opposed to any active computation
on the encoded data). We think of the code during the error part of the
protocol as sitting idle in the presence of an undesirable environment. The
effect of the environment is what produces errors in the code, so the channel
$\mathcal{E}$ implicitly contains information about the nature of the
environment and how long the code is exposed to the environment. We denote
such an error channel by $\mathcal{E}_{t}$, where $t$ corresponds to the time
for which the code is exposed to the environment. Correspondingly, the
recovery operation $\mathcal{R}$ may depend on $t$ and indeed on the whole
structure of $\mathcal{E}_{t}$.

In this paper we consider the environment to be a thermal reservoir at some
temperature $T$. There is still freedom in determining exactly how the
reservoir couples to the system, but from an effective theory perspective we
should expect all allowed couplings to be present with some strength. In this
context, the need for active error correction can be phrased like this. Given
a code $\mathcal{C}$, the code requires active error correction if, as
$t\rightarrow\infty$, there is no $\mathcal{R}$ that can correct
$\mathcal{E}_{t}$. In other words, a code requires active error correction if,
given enough time, the error channel will eventually degrade the information
beyond the ability of any physically conceivable recovery operation to repair.

More generally, we can define a memory time $t_{\text{mem}}$ as the time
beyond which $\mathcal{E}_{t}$ degrades the information beyond any hope of
recovery. Of course, this time depends on the code and the error model, but it
also depends on what recovery operations we can conceivably perform and on how
well we want to be able to recover the information (e.g. the degree of
approximation in our approximate code). A perfect passive memory is one for
which $t_{\text{mem}}=\infty$. Typically, we can only hope for this sort of
behavior in a thermodynamic limit, so we should really view a passive memory
as a sequence of codes indexed by the physical size $n$ for which
$\lim_{n\rightarrow\infty}t_{\text{mem}}(n)=\infty$. Even more practically,
what really matters is how $t_{\text{mem}}$ scales with $n$ and other
parameters of the code.

For the purpose of potentially implementing a quantum memory in the lab, we
may also want to impose some additional rules of the game. One variant of such
rules are called the Caltech rules [20, 21]. For example, we might obtain a
long memory time by sending the energy cost of errors to infinity “by hand”,
but this could be considered “cheating” unless there is a good physical reason
for this large energy penalty. The idea is that we indeed want to effectively
forbid errors, but this should happen naturally without sending any
microscopic energy scale to infinity. Of course, if one does have a system
where errors are particularly costly, then it might form a good enough passive
memory in the sense that $t_{\text{mem}}$ is long enough for practical
purposes. To give another example, to be able to construct codes in physical
three-dimensional space, we might want some restrictions on the kinds of codes
we consider, e.g. trying to avoid the need for too many long-range
interactions [21]. We note that all the models we consider here do have non-
local interactions between the large $N$ degrees of freedom333Note that it is
somewhat ambiguous what constitutes the local degrees of freedom in the large
$N$ gauged oscillator model as they are not simply tensor product of qubits.
We show in later sections that there is a construction where two such gauged
oscillators are spatially local with respect to each other. However, each one
contains a large $N$ degrees of freedom even though there is no sense in which
they occupy an extended region in space. It is in this sense that these
degrees of freedom are related to each other non-locally..

###  1.2 Background on quantum gauge theories

In this subsection, we review the key concepts of quantum gauge theories
needed for our discussion. This section is aimed at readers with a more
quantum information background.

At the most basic level, a gauge theory is defined by giving a gauge group
$\cal{G}$ and a set of degrees of freedom transforming under $\cal{G}$ with
the requirement that all physical states and operators are invariant under
$\cal{G}$. This invariance condition is known as Gauss’ law, in honor of its
roots in electromagnetism. An immediate question is why one would use such a
redundant description at all? Experience shows that one can always use a
gauge-invariant description of the states and operators, but such a
description often hides important underlying features of the physics. For
example, loop-like operators known as Wilson lines provide a class of gauge
invariant observables in quantum field theories with gauge fields, but the
extended nature of these operators can obscure aspects of locality in the
theory. In special cases, one can even have an equivalence between a local
gauge theory and another local theory with no explicit gauge structure, for
example, the duality between the 3d statistical Ising model and 3d $Z_{2}$
gauge theory.

We study such gauge theories in large part because they form a crucial part of
our description of nature. The quantum theory of electromagnetism is a gauge
theory where the gauge group is related to the abelian Lie group $U(1)$. The
more general Standard Model is also a gauge theory with a non-Abelian gauge
group of the form $SU(3)\times SU(2)\times U(1)$. Gauge theories also arise in
condensed matter systems, for example, when describing so-called spin liquid
phases, and, as we mentioned above, they also have interesting error
correction properties.

Apart from the motivations from known models of high-energy particle physics
and condensed matter systems, another motivation to consider gauge theories
comes from quantum gravity. Since the discovery of AdS/CFT, it is generally
believed that certain $d$ dimensional $SU(N)$ gauge theories with large $N$
admit a dual geometric description in terms of a $d+1$ dimensional theory of
gravity. From a quantum information perspective, this correspondence provides
the possibility of hiding quantum information deep inside the emergent
geometry. This has inspired a lot of research in the field of holographic
quantum error correcting codes [7], with explicit models such as the HaPPY
code [11], the effective field theory code [22], and many other related
constructions using tensor networks [23, 24, 25]. Current holographic codes
often utilize the extra dimension, such that any local/UV errors can not
access the information deep inside along the extra radial direction. However,
very little of the error correcting power in these constructions come directly
from the gauge symmetries in AdS/CFT. In this work, we will focus on a
different kind of protection from the large $N$ gauge symmetries, which can be
related to the small Newton’s constant in these holographic theories [16, 10].
This opens another possibility to design a new kind of quantum memory.

An interesting structure we exploit in the large $N$ gauge theory is its
sparse density of excitations at low energy. Due to Gauss’ law, the only
allowed excitations of a $SU(N)$ gauge theory on any compact manifold are
gluons that combine into gauge singlets. Therefore, the density of states is
determined by counting the number of gauge singlets in a given energy level.
It is a general feature that the density of states in such model scales
exponentially with some power of energy and does not depend on $N$. As is
pointed out in [26], even for a free gauge theory with large $N$, the
projection to the gauge invariant subspace induces effective interactions
among the gluons, which is negligible at high energy but dominates at the low
energy regime. In fact, for a weakly coupled $SU(N)$ gauge theory, there is a
phase transition analogous to the comfinement/deconfinement transition in QCD,
in the limit of $N\rightarrow\infty$. As we show, this induces a critical
temperature, $T_{c}$, below which the quantum information can be long-lived
but above which the information is quickly destroyed. In the geometric
description, this is believed to be dual to the Hawking-Page transition and
the high temperature phase can be thought of as a black hole that kills the
quantum information [10].

To zero in on just the effects of gauge symmetry, we distill only the
necessary ingredients from the above gauge theory intuitions and consider a
gauge theory with no spatial geometry at all. It turns out that even the
simplest $0+1$-dimensional gauge theory of harmonic oscillators is sufficient
for our purposes. The inverted harmonic oscillator has been studied as a dual
to a two dimensional string theory and can be considered as toy model of
holography [27, 28, 29, 30]. Here we start from a similar but non-inverted
harmonic oscillator model and try to build a quantum memory from it.

###  1.3 Organization and overview of results

The rest of this paper is organized as follows. In Section 2, we describe in
detail the large $N$ matrix quantum mechanics model with $SU(N)$ gauge
symmetry. We then discuss how it may be construed as an approximate error
correcting code and define various bath models we will use for the rest of the
paper. We then discuss the mutual information diagnostic for memory time and
distill the core conclusions of our work in the form of Theorem 1. It shows
that if a large $N$ system has a sparse low energy spectrum, couples to the
bath uniformly, and is an approximate QECC, then it has a memory time
polynomial in $N$ at low enough temperatures. In the ensuing Sections, we
deploy this theorem to various quantum memory constructions.

In Section 3, we study two concrete examples of quantum memories built from
two gauged harmonic oscillator models coupled to a thermal reservoir. The
memory time scales as $N^{2}$ in both models. For the first model, the $N^{2}$
life time persists even at high temperature, but the model requires non-local
logical operations. We then construct a local model which has a $N^{2}$ memory
time in the low temperature phase, but requires a large coupling of strength
$\sim\log N$ between the gauged oscillators. At zero temperature limit
$T\rightarrow 0$, the first model has a constant $N^{2}$ life time, while the
second model has a divergent life time that scales as
$N^{2}e^{\frac{\omega}{T}}$, where $\omega$ is some energy scale of this
model.

In Section 4, we consider what happens if the gauge constraint is not imposed
exactly, but only as an energetic cost in the Hamiltonian. A large number of
non-singlet excitations now modify the physics, but the information retention
time can be still be quadratic in $N$ if the energy scale of non-singlet
excitations is large.

Finally, in Section 5, we ask to what extent does robustness depend on the
symmetry being exact and the Hilbert space of quantum systems being infinite
dimensional. On the practical level, we want the quantum memories to be
implementable on finite dimensional quantum systems like qubits or qudits. To
this end, we consider an even more relaxed situation where the $SU(N)$
symmetry itself is not only ungauged but also approximate. We then build a
quantum memory consisting of $N^{2}$ number of physical qudits that interact
non-locally with each other. Qualitatively, this model behaves similarly to
the one in Section 4 in the low energy sector, except small correction terms
that suppressed to order at most $\sim\log N/N$. We find that if we require
the total energy of the interaction scales as $\log{N}$ and a reasonable
pattern of coupling strengths, then Theorem 1 again holds to guarantee a
$N^{2}$ memory time. The zero temperature limit of life times for the last two
models also diverge as $N^{2}e^{\frac{\omega}{T}}$.

##  2 Matrix model and set up

###  2.1 Matrix Quantum Mechanics

We wish to study a gauged matrix model of harmonic oscillators. First, recall
the classical Lagrangian of a single harmonic oscillator,

$L=\frac{1}{2}(\dot{x}^{2}-\omega^{2}x^{2}).$ (3)

The canonical momentum is $p=\frac{\partial L}{\partial\dot{x}}=\dot{x}$, in
terms of which we also have the Hamiltonian description,

$H=\frac{1}{2}(p^{2}+\omega^{2}x^{2}).$ (4)

The system is quantized in the usual way by setting $[x,p]=i\hbar$.

A generalization can be made by promoting the variable $x(t)$ to an $N\times
N$ matrix $X(t)$ and promoting products of variables to matrix multiplication.
We hasten to emphasize that this $N\times N$ matrix structure does not refer
to a quantum Hilbert space; rather it is a way of organizing $N^{2}$ variables
into a single unit. Then one gets an action for a matrix harmonic oscillator,

$L=\frac{1}{2}\operatorname{Tr}[(\partial_{t}X)^{2}-\omega^{2}X^{2}].$ (5)

Again, let us emphasize that the trace is not a quantum Hilbert space trace
but rather a sum over the matrix degrees of freedom. The diagonal elements of
$X$ are real variables, while the off-diagonal components are complex
variables coming in conjugate pairs.

The momentum conjugate to $X^{i}_{j}$ is $P^{j}_{i}=\dot{X}^{j}_{i}$. We
quantize the theory by specifying canonical commutators for all pairs of
conjugate variables. If we extend quantum Hermitian conjugation to include a
tranpose of the matrix degrees of freedom, then the quantum variable $X$ is
Hermitian. In terms of its matrix elements $X^{i}_{j}$, Hermiticity implies
that $X^{i}_{i}$ are quantum Hermitian while $X^{i}_{j}$ and $X^{j}_{i}$ are
quantum Hermitian conjugates, $(X^{i}_{j})^{\dagger}=X^{j}_{i}$. Note that
$X^{i}_{j}$ is not quantum Hermitian when $i\neq j$, but the above formula is
still correct.

After quantization, it is convenient to work in terms of creation and
annihilation operators,

$\begin{split}a^{i}_{j}&=\sqrt{\frac{\omega}{2}}(X^{i}_{j}+iP^{i}_{j})\\\
a^{\dagger j}_{i}&=\frac{1}{\sqrt{2\omega}}(X^{j}_{i}-iP^{j}_{i}),\end{split}$
(6)

which satisfy

$[a^{i}_{j},a^{\dagger k}_{l}]=\delta^{i}_{l}\delta^{k}_{j}.$ (7)

The Hamiltonian written in terms of $a$ and $a^{\dagger}$’s takes the simple
form

$H=\omega\operatorname{Tr}(a^{\dagger}a)$ (8)

after subtracting the zero-point energy. There is a unique vacuum state
$|0\rangle$ satisfying

$a^{i}_{j}|0\rangle=0,~{}\forall i,j.$ (9)

Acting creation operators on this vacumm then generates all the Fock states
that span the Hilbert space.

Looking back at the Lagrangian in Eq. (5), we see that it has a symmetry under
unitary transformations $X(t)\rightarrow UX(t)U^{\dagger}$, with $U\in U(N)$.
(This is why we consider complex off-diagonal entries in $X$.) This can be
made into a local symmetry $U(t)$, meaning a symmetry where $U(t)$ can be
different at each different time, as follows. We add a matrix gauge field
$A(t)$ and promote the time derivative to a covariant derivative, defined as

$D_{t}X=\partial_{t}X-[A,X]$ (10)

with $A(t)$ transforming as $A(t)\rightarrow
U(t)A(t)U(t)^{\dagger}-iU(t)\partial_{t}U(t)^{\dagger}$. One can check that
the covariant derivative transforms in the same way as $X(t)$, which is
$D_{t}X\rightarrow U(t)D_{t}XU(t)^{\dagger}$. Therefore, the Lagrangian is
invariant under this transformation. We call this symmetry with respect to the
time-dependent transformation the gauge symmetry.

The model equipped with such a symmetry can be regarded as a redundant
description of a physical system whose Hilbert space is defined to include
only gauge-invariant states. Indeed, in this model, $A(t)$ can be view as a
Lagrange multiplier. Fixing a gauge $A=0$, and integrating out $A$ in the path
integral representation gives the Gauss law constraint:

$G=[X,\dot{X}]=0.$ (11)

This condition is implemented by requiring that physical states $|\psi\rangle$
be annihilated by $G$, $G|\psi\rangle=0$.

One can see that the states satisfying this constraint are the symmetry
singlets, generated by traces of operators acting on the vacuum,

$|{\psi}\rangle=\operatorname{Tr}(a^{\dagger
n_{1}})\operatorname{Tr}(a^{\dagger n_{2}})\cdots\operatorname{Tr}(a^{\dagger
n_{k}})|{0}\rangle.$ (12)

These states, after being normalized, span a subspace of the ungauged matrix
model’s Hilbert space. They are also eigenstates of the Hamiltonian (8) such
that

$H|\psi\rangle=\omega\sum_{i=1}^{k}n_{i}|{\psi}\rangle.$ (13)

###  2.2 Physical states and the Knill-Laflamme Condition

We now describe in broad terms how the properties of these gauge-singlet
states can be leveraged to generate approximate codes obeying an approximate
version of the KL condition. First some terminology. If a gauge invariant
operator contains only a single trace, it is called single trace operator.
Operators with multiple traces are called multi-trace. Using the commutation
relation, Eq. (7), one can derive a useful commutator rule for single trace
operators, with $n_{1,2}$, $m_{1,2}$ being $O(1)$ number,

$[\operatorname{Tr}(a^{n_{1}}a^{\dagger
m_{1}}),\operatorname{Tr}(a^{n_{2}}a^{\dagger
m_{2}})]=(n_{1}m_{2}-n_{2}m_{1})\operatorname{Tr}(a^{n_{1}+n_{2}-1}a^{\dagger
m_{1}+m_{2}-1})+O\left(\frac{1}{N^{2}}\right).$ (14)

Using these rules, overlaps of states generated by acting multi-trace
operators on the ground state can be obtained to leading order in large $N$.
For example, the squared norm of the state
$\operatorname{Tr}(a^{\dagger})|0\rangle$ is

$\langle 0|\operatorname{Tr}(a^{n})\operatorname{Tr}(a^{\dagger
n})|0\rangle=\langle 0|[\operatorname{Tr}(a^{n}),\operatorname{Tr}(a^{\dagger
n})]|0\rangle=nN^{n}+O\left(\frac{1}{N^{2}}\right).$ (15)

We also find that states with the same energy are almost orthogonal, with
overlaps suppressed by a power of $N$. For example,

$\langle
0|\frac{\operatorname{Tr}(a^{m})\operatorname{Tr}(a^{n-m})}{N^{\frac{m}{2}}}\frac{\operatorname{Tr}(a^{\dagger
n})}{N^{\frac{n}{2}}}|0\rangle=O\left(\frac{1}{N}\right).$ (16)

Moreover, one can show that a factorization property holds for products of
single trace operators. As an example, given four single trace operators with
vanishing expectation value, $O_{i}$, $O_{j}$, $O_{k}$ and $O_{l}$, we have

$\begin{split}\langle
0|O_{i}O_{k}O_{l}^{\dagger}O_{j}^{\dagger}|0\rangle&\sim\langle
0|O_{i}O_{j}^{\dagger}|0\rangle\langle 0|O_{k}O_{l}^{\dagger}|0\rangle+\langle
0|O_{i}O_{l}^{\dagger}|0\rangle\langle
0|O_{k}O_{j}^{\dagger}|0\rangle+O\left(\frac{1}{N^{2}}\right)\\\ \langle
0|O_{i}O_{k}O_{j}^{\dagger}|0\rangle&\sim
O\left(\frac{1}{N}\right).\end{split}$ (17)

If we could find a way to classify $O_{i}$ and $O_{j}$ as logical operators
and $O_{k}$ and $O_{l}$ as errors, then this equation would be interpreted as
an approximate version of Knill-Laflamme condition (aKL). While one may
deviate from KL in many forms, the approximate KL condition we consider in
this work could be written as

$\begin{split}\langle
i|E_{a}^{\dagger}E_{b}|j\rangle&=f_{ab}\delta_{ij}+\frac{g^{ij}_{ab}}{N^{2}}\\\
\langle i|E_{a}|j\rangle&=\frac{e_{a}^{ij}}{N}.\end{split}$ (18)

In the above formulae, we are imagining that $E_{a}$ and $E_{b}$ are some
general non-identity error operators. In general, we expect them to be any
gauge invariant operators that satisfy some constraints, depending on how we
couple the model to a thermal environment. We will give concrete examples in
Section 3, and explicitly calculate the functions $f$, $g$ and $e$. These
examples will show such large $N$ matrix quantum mechanics models can used to
construct approximate quantum error correcting codes [14].

For later use, we also briefly discuss the counting of gauge invariant
operators. If we have only one matrix harmonic oscillator, the general
invariant states are in the form of Eq. (12). At energy level $n\omega$, the
number of distinct states (which are also approximately orthogonal) is equal
to the number of partitions of integer $p(n)$. For example, at level $n=3$.
The three partitions are,

$\begin{split}&3=3+0\\\ &3=2+1\\\ &3=1+1+1.\end{split}$ (19)

Correspondingly, we have three gauge invariant states at energy $3\omega$,
which are

$\begin{split}\operatorname{Tr}(a^{\dagger 3})|0\rangle,\ \
\operatorname{Tr}(a^{\dagger 2})\operatorname{Tr}(a^{\dagger})|0\rangle,\ \
\operatorname{Tr}(a^{\dagger})^{3}|0\rangle.\end{split}$ (20)

In general, we can have arbitrary $k\leq n$ number of matrix harmonic
oscillators, each carrying energy that is an integer multiple of $\omega$. We
label the individual matrix operator by its subindex, $a_{1}$, $a_{2}$,
$a_{3}$, $\cdots$. The counting becomes more complicated with multiple
oscillators, so let us focus on the states generated by single trace
operators. At energy level $n\omega$, these states can be written as

$\operatorname{Tr}(P\\{a_{1}^{\dagger n_{1}}a_{2}^{\dagger n_{2}}\cdots
a_{k}^{\dagger n_{k}}\\})|0\rangle,$ (21)

with $n=\sum_{i}n_{i}$. The state $|0\rangle$ is tensor product of $k$ ground
states of each oscillator. We use the symbol $P$ to denote a particular
permutation of all the operators inside the bracket. By counting how many
different ways of splitting $n$ and the number of permutations, one can show
that the number of these single trace operators asymptotically scales as
$k^{n}$ when $n$ is large enough. For simplicity, in the concrete models that
will be discussed in the following sections, we consider no more than two
matrix modes.

###  2.3 Bath model

We have seen that if errors and logical operators can be chosen properly, then
an approximate Knill-Laflamme condition holds. Based on this, we expect that a
single error operator can only corrupt the logical information at order
$\frac{1}{N^{2}}$. At the same time, gauge invariance constrains the number of
possible errors that may occur in the form of excitations. Since the density
of states for the gauged oscillator grows with energy but is independent of
$N$ below some critical value, we may expect such a system to serve as a good
quantum memory even at finite temperature. To show this explicitly, we couple
the proposed quantum memory to a thermal reservoir and calculate the memory
time $t_{\text{mem}}$.

We make the standard assumption that the bath is Markovian and couples to the
system locally. More specifically, we model the bath as collection of bosonic
modes $b_{l}$, e.g. a collection of oscillators, with the following thermal
spectral function [31],

$\operatorname{Tr}(\rho_{B}b_{k}(\nu)b_{l}^{\dagger}(-\nu))=\delta_{kl}\left|\frac{\nu}{1-e^{-\beta\nu}}\right|:=\delta_{kl}\gamma(-\nu).$
(22)

Physically, this means that the probability of having a given error that adds
energy $\epsilon=-\nu$ to the system is suppressed by a thermal factor
$\gamma(\epsilon)$.

We will discuss three scenarios for the structure of the bath. In particular,
there is a degree of ambiguity in what constitutes a generic physical bath for
these systems, so we consider a number of models of varying levels of
generality.

Case 1: Consider two matrix harmonic oscillators, $a_{1}$ and $a_{2}$,
separated by a large distance in space. The reason for introducing two such
oscillators will become clear in the next section. We model this setup by
assuming that these oscillators couple to the thermal bath independently with
order one coupling constants. This means we assume that thermal errors can
include any multi-trace operators provide they are uncorrelated between
$a_{1}$ and $a_{2}$. By uncorrelated error, we mean that the modes $a_{1}$ and
$a_{2}$ can not appear inside the same trace; instead errors must be of the
form

$E_{L,R}=\prod_{(r_{1},r_{2})\in L,\ (s_{1},s_{2})\in
R}\frac{:\operatorname{Tr}(P\\{a_{1}^{\dagger
r_{1}}a_{1}^{r_{2}}\\})::\operatorname{Tr}(P\\{a_{2}^{\dagger
s_{1}}a_{2}^{s_{2}}\\}):}{N^{\frac{r_{1}+r_{2}+s_{1}+s_{2}}{2}}}$ (23)

with $L$ and $R$ being some sets of integer 2-tuples. The symbol $P$ denotes
the sum over all distinct arrangements of the creation and annihilation
operators inside the bracket. The symbol $::$ means normal ordering and the
$N$-dependent normalization is to ensure that the single trace operator of
each $a$, $a^{\dagger}$ arrangement has order one norm (see Eq. (15)). In
terms of these errors, the thermal Hamiltonian is

$H_{\rm thermal}=\sum_{L,R}\lambda_{L,R}b_{L,R}E_{L,R}+h.c.,$ (24)

where h.c. stands for hermitian conjugate. We will discuss this case in more
detail in Section 3.1.

Case 2: We assume the coupling with the bath modes $b_{l}$ are through single-
trace operators, with order one coupling constants $\lambda_{l}$’s. We still
consider modes $a_{1}$ and $a_{2}$ from two different matrix oscillators. The
thermal coupling Hamiltonian is in the form of

$H_{\rm
thermal}=\sum_{\\{n_{k}\\}}\lambda_{\\{n_{k}\\}}b_{\\{n_{k}\\}}\frac{:\operatorname{Tr}(P\\{a_{1}^{\dagger
n_{1}}a_{2}^{\dagger
n_{2}}a_{1}^{n_{3}}a_{2}^{n_{4}}\\}):}{N^{\frac{\sum_{k=1}^{4}n_{k}}{2}}}+h.c.,$
(25)

At first sight, it appears that coupling only to single-trace operators is a
strong assumption, because in general we should allow all gauge invariant
operators to couple with the bath. However, this is still general enough since
single-trace errors accumulate into multi-trace errors as time increases. We
will discuss this scenario in more detail in Section 3.2.

Case 3: We consider the most general coupling to the thermal bath where each
of the $N^{2}$ harmonic oscillator modes can couple to the bath independently
without requiring the coupling terms be gauge invariant. In addition, we
assume that the coupling constant decays exponentially with the number of
creation and annihilation operators involved. This is a reasonable assumption
because the correlated errors that involves more operators occur with smaller
probability. We explain this in detail in Section 4, where we consider the
ungauged model.

###  2.4 Mutual information diagnostic

Given the assumptions on the bath described above, the time evolution of our
system can be described by a Lindblad equation after tracing out the bath,

$\dot{\rho}(t)=-i[H,\rho(t)]+\sum_{a}\lambda_{a}(\nu_{a})E_{a}\rho(t)E_{a}^{\dagger}-\frac{1}{2}\sum_{a}\lambda_{a}(\nu_{a})\\{E_{a}^{\dagger}E_{a},\rho(t)\\},$
(26)

where $\lambda_{a}(\nu_{a})=\lambda_{a}\gamma(\nu_{a})$, is the coupling
constant multiplied with thermal factor. $\nu_{a}$ is the energy of
excitation, and the $\\{E_{a}\\}$ depend on the interaction Hamiltonian.

The generic expectation is that the above system dynamics will lead to rapid
information loss within the system. We wish to show that this is not the case
for the right kind of code constructed from the matrix degrees of freedom. To
quantitatively measure the memory time, we entangle the encoded
$d$-dimensional system with a reference $R$. The initial state is thus

$\begin{split}&\frac{1}{\sqrt{d}}\sum_{i=1}^{d}|i\tilde{i}\rangle:=\frac{1}{\sqrt{d}}\sum_{i=1}^{d}|i\rangle_{S}|i\rangle_{R}\\\
&\rho_{0}=\frac{1}{d}\sum_{i,j=1}^{d}|i\tilde{i}\rangle\langle
j\tilde{j}|,\end{split}$ (27)

and the correlation between the reference is maximal.

The goal is then to understand how the system-reference correlation evolves
with time. The state evolves under the Lindblad dynamics into

$\begin{split}\rho(t)&=e^{\mathcal{L}t}\rho_{0},\\\
\textrm{where}~{}~{}\mathcal{L}\mathcal{O}&=\sum_{k}\lambda_{k}E_{k}\mathcal{O}E_{k}^{\dagger}-\sum_{k}\frac{\lambda_{k}}{2}\\{E_{k}^{\dagger}E_{k},\mathcal{O}\\}.\end{split}$
(28)

In terms of the system and reference reduced density matrices
$\rho_{S}(t)=\operatorname{Tr}_{R}[\rho(t)]$ and
$\rho_{R}(t)=\operatorname{Tr}_{S}[\rho(t)]$, we may quantify the correlation
using the mutual information,

$I(S:R)(t)=S(\rho_{S}(t))+S(\rho_{R}(t))-S(\rho(t)).$ (29)

When the mutual information is maximal, $I=2\ln d$, we say that the quantum
information is stored in the quantum memory. However, as $I(S:R)$ decays, the
information of the encoded qubit leaves the memory and dissipates into the
environment. For the sake of simplicity, we will focus on quantum memory
constructions that store a single encoded qubit in the ensuing sections, so
$d=2$.

We now tie together all the ingredients above in a theorem. In general,
different couplings with the bath can induce different errors on the system.
However, as long as such errors are not too numerous and can be suppressed by
the code, i.e., satisfying the aKL condition, then the following theorem
guarantees the memory time of the quantum information is relatively long. To
be precise, this theorem shows that the mutual information diagnostic decays
polynomially with $t$ at a rate suppressed by $\frac{1}{N^{2}}$, which implies
that the memory time can scale polynomially with $N$.

######  Theorem 1.

Consider a system whose low energy states form an approximate quantum error
correcting code that encodes a logical qudit $S$ and a set of physical errors
$\mathcal{E}$ induced by the system-bath coupling (c.f. Section 2.3). Provided
the following conditions are satisfied,

  1. 1.

(Sparse spectrum) There exists an energy $\epsilon_{0}$ below which the
effective number of errors444These are errors that contribute non-trivially in
the Lindblad equation. in $\mathcal{E}$ with energy $\epsilon\leq\epsilon_{0}$
is independent of $N$ and bounded by $\exp(\mu\epsilon)$ for some $\mu>0$,

  2. 2.

(Uniform coupling) The thermal Hamiltonian that couples each $E\in\mathcal{E}$
to an independent bath operator has coupling constant no larger than $O(1)$,

  3. 3.

(Approximate error correction) $\forall E_{a},E_{b}\in\mathcal{E}$, the
approximate Knill-Laflamme condition (Eq. (18)) is satisfied,

then the mutual information $I(S,R)$ defined in Eq. (29) is given by

$I(S:R)(t)=2\ln{d}-K\left(\frac{t}{N^{2}}\right),$ (30)

where $K(x)$ is a polynomial function of $t$, for temperature
$T<\frac{1}{\mu}$.

###### Proof.

The proof is given in Appendix A. ∎

In the following sections, we will construct different systems coupled via
different bath models, as outlined in Section 2.3, and prove that they satisfy
the assumptions in the theorem, such as the approximate Knill-Laflamme
condition (aKL) for a corresponding set of error operators and logical states.

Remark:  Note that the above theorem does not require the code subspace to be
the ground space of the system. Therefore, the $N^{2}$ scaling of memory time
holds generally for these models no matter how the code subspace is chosen as
long as the requisite conditions are satisfied. Intuitively, this is because
the transition rate between states are $1/N^{2}$ suppressed. However, as
excited states can still slowly decay into lower energy states at the same
rate, their memory times need not be longer by lowering the temperature.

On the other hand, for quantum memories where we also choose the code subspace
to be the ground space, the memory lifetime can be further extended by a
multiplicative factor $\exp(\omega/T)$, which is related to the effect of
thermal suppression. This renders the total low temperature memory time to be
of $O(N^{2}\exp(\omega/T))$ for some energy scale $\omega$ as we see in some
of the examples below.

##  3 Fully gauged model

###  3.1 Non-local model

As the first example, consider two gauged harmonic oscillators $a_{1}$ and
$a_{2}$ with $SU(N)$ gauge symmetry, separated by a large distance. Their
interaction is weak, and each couples to the thermal reservoir locally and
independently555More realistically, thermal noise can have finite correlation.
In such cases, we assume that the distance separating the two oscillators is
much larger than the correlation length so our assumptions remain good
approximations.. The Hamiltonian takes the following form,

$H=\omega\operatorname{Tr}(a_{1}^{\dagger}a_{1})+\omega\operatorname{Tr}(a_{2}^{\dagger}a_{2})+\sum_{L,R}\lambda_{L,R}b_{L,R}E_{L,R}+h.c.,$
(31)

where $E_{L,R}$’s are the possible errors acting on this model, and
$\lambda_{L,R}$’s are $O(1)$ coupling constants.

The general gauge invariant states are given in Eq. (21) with $k=2$. However,
since we assume that these modes only couple to the thermal bath locally, the
set of possible error operators is a proper subset of all physical operators.
These errors are uncorrelated in the sense that each error operator contains
only a single trace either the mode $a_{1}$ or the mode $a_{2}$, but never
both. Given two sets of integer 2-tuples $L$ and $R$, we can define the errors
as,

$E_{L,R}=:\prod_{(r_{1},r_{2})\in L,\ (s_{1},s_{2})\in
R}\frac{\operatorname{Tr}(P\\{a_{1}^{\dagger
r_{1}}a_{1}^{r_{2}}\\})\operatorname{Tr}(P\\{a_{2}^{\dagger
s_{1}}a_{2}^{s_{2}}\\})}{N^{\frac{r_{1}+r_{2}+s_{1}+s_{2}}{2}}}:.$ (32)

where the symbol $P$ denotes the sum over all distinct arrangements of the
$a_{i},a_{i}^{\dagger}$ operators666For example, $a^{\dagger 2}a^{2}$ and
$a^{\dagger}aa^{\dagger}a$ are counted as two distinct arrangements. . In
general, $P\\{a^{\dagger r_{1}}a^{r_{2}}\\}$ admits
$(r_{1}+r_{2})!/r_{1}!r_{2}!$ terms of distinct arrangements. Each operator is
brought to normal ordering, and divided by some power of $N$ to make sure that
these operators have $O(1)$ norm.

#### Logical states

Given this set of error operators, we can construct logical states that are
protected against them. In general, we may choose the logical states to be

$|\mathcal{I}\rangle:=\prod_{(m,n)\in\mathcal{I}}\frac{\operatorname{Tr}(a_{1}^{\dagger
m}a_{2}^{\dagger n})}{N^{\frac{m+n}{2}}}|0\rangle_{12}.$ (33)

with $m,n\geq 1$ and belonging to some set $\mathcal{I}$ of integer 2-tuples.
The state $|0\rangle_{12}$ is the tensor product of the vacuum states of the
two oscillators. For simplicity, let us focus on a code subspace spanned by
two states and treat it as a logical qubit777The logical subspace we use is
not the lowest energy subspace of the Hamiltonian, so they have the tendency
to decay into lower energy states. However, the decay rate is suppressed by
$\frac{1}{N^{2}}$, as discussed in Appendix B.1.:

$\begin{split}|\tilde{1}\rangle&=\frac{\operatorname{Tr}(a_{1}^{\dagger
2}a_{2}^{\dagger 2})}{N^{2}}|0\rangle_{12}\\\
|\tilde{2}\rangle&=\frac{\operatorname{Tr}(a_{1}^{\dagger}a_{2}^{\dagger})^{2}}{\sqrt{2}N^{2}}|0\rangle_{12}.\end{split}$
(34)

They are not orthogonal states, since
$\langle\tilde{1}|\tilde{2}\rangle=\frac{\sqrt{2}}{N}+O(\frac{1}{N^{3}})$.
However, it is straightforward to identify an orthonormal basis
$|\tilde{\uparrow}\rangle$ and $|\tilde{\downarrow}\rangle$, defined as

$\begin{split}&|\tilde{\uparrow}\rangle=|\tilde{1}\rangle\\\
&|\tilde{\downarrow}\rangle=\frac{|\tilde{2}\rangle-\langle\tilde{2}|\tilde{1}\rangle|\tilde{1}\rangle}{\sqrt{1-\langle\tilde{2}|\tilde{1}\rangle^{2}}}.\end{split}$
(35)

#### Errors and aKL

Because we assume errors are gauge invariant, they are expressible as single
or multi-trace operators (see Eq. (32)). For the purpose of invoking Theorem
1, it suffices to check that aKL holds for these errors,

$\begin{split}\langle\tilde{i}|E_{L^{\prime},R^{\prime}}^{\dagger}E_{L,R}|\tilde{j}\rangle&=f_{L^{\prime}R^{\prime},LR}\delta_{ij}+\frac{g_{L^{\prime}R^{\prime},LR}^{ij}}{N^{2}}+O(\frac{1}{N^{4}})\\\
\langle\tilde{i}|E_{L,R}|\tilde{j}\rangle&=\frac{e_{L^{\prime}R^{\prime},LR}^{ij}}{N}+O(\frac{1}{N^{3}}).\end{split}$
(36)

To build up some intuition, let us first examine a few instructive examples
before we arrive at the more general relation (36). To start, one can check
that the operator $\frac{\operatorname{Tr}(a_{1}^{\dagger}a_{1})}{N}$ has the
logical states as eigenstates,

$\frac{\operatorname{Tr}(a_{1}^{\dagger}a_{1})}{N}|\tilde{i}\rangle=\frac{2}{N}|\tilde{i}\rangle.$
(37)

Therefore, it is a logical identity operator. We can similarly check that the
Hamiltonian that generates the system dynamics is proportional to the logical
identity, which is necessary for a memory.

Next, let us examine a class of single-trace errors of the form
$E_{n}=\frac{\operatorname{Tr}(a_{1}^{\dagger n})}{\sqrt{n}N^{\frac{n}{2}}}$.
As a potential bit flip error, we check the transition amplitude,

$\langle\tilde{\downarrow}|E_{m}^{\dagger}E_{n}|\tilde{\uparrow}\rangle=O(\frac{1}{N^{3}}).$
(38)

The bit flip error only carries $n$ dependence at order $O(\frac{1}{N^{3}})$.
Similarly, the phase errors are also suppressed,

$\begin{split}&\langle\tilde{\uparrow}|E_{m}^{\dagger}E_{n}|\tilde{\uparrow}\rangle=\delta_{nm}\left[1+\frac{3n}{N^{2}}+O\left(\frac{1}{N^{4}}\right)\right]\\\
&\langle\tilde{\downarrow}|E_{m}^{\dagger}E_{n}|\tilde{\downarrow}\rangle=\delta_{nm}\left[1+\frac{2n}{N^{2}}+O\left(\frac{1}{N^{4}}\right)\right].\end{split}$
(39)

Note that the leading order terms of these equations come from factorized two
point functions $\langle
0|E^{\dagger}_{m}E_{n}|0\rangle\langle\tilde{i}|\tilde{j}\rangle$, and they
represent the non-connected contribution to the four point correlator. The
sub-leading corrections are from the connected part, shown diagrammatically in
Figure 1. It is this part that induces phase and bit-flip errors between the
two logical states. In a similar way, one can check that the aKL condition
holds for other single-trace errors. We present the more comprehensive
calculations in Appendix B.1.

$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger
2})$$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger
2})$$E_{n}$$E_{n}$$+$$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger
2})$$E_{n}$$E_{n}$$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger 2})$
Figure 1: The inner product can be viewed as scattering between two particles.
$(1,0)$ represents the logical state $|\tilde{\uparrow}\rangle$. The connected
part is suppressed by $\frac{1}{N^{2}}$.

The conclusion of this analysis is that single-trace operators can only
produce $1/N$ suppressed errors, thus giving rise to the approximate Knill-
Laflamme condition (aKL). Using the results for single-trace errors, one can
also derive the aKL conditions for any multi-trace errors. As an example,

$\begin{split}&\langle\tilde{\uparrow}|E_{p}^{\dagger}E_{q}^{\dagger}E_{m}E_{n}|\tilde{\uparrow}\rangle=F(m,n,p,q)+\frac{G_{\uparrow}(m,n,p,q)}{N^{2}}+O(\frac{1}{N^{4}})\\\
&\langle\tilde{\downarrow}|E_{p}^{\dagger}E_{q}^{\dagger}E_{m}E_{n}|\tilde{\downarrow}\rangle=F(m,n,p,q)+\frac{G_{\downarrow}(m,n,p,q)}{N^{2}}+O(\frac{1}{N^{4}}).\end{split}$
(40)

The function $F$ is the same in both correlators, so it represents a
correctable error. In contrast, different $G$ functions appear in the two
equations. Their difference, $G_{\uparrow}-G_{\downarrow}$, represents a small
uncorrectable error. They both have simple diagrammatic representation as
shown in Figure 2 and Figure 3, from which we can read off the function
$F(m,n,p,q)$. It must take the form

$\begin{split}&F(m,n,p,q)=(\delta_{m}^{p}\delta_{n}^{q}+\delta_{m}^{q}\delta_{n}^{p})+\delta_{m+n}^{p+q}\frac{f(p,q,m,n)}{N^{2}},\end{split}$
(41)

where the $O(\frac{1}{N^{2}})$ term comes from the partly connected diagram in
Figure 2. Since $F$ does not depend on the logical information and will not
contribute to uncorrectable errors, it is not pertinent to the subsequent
discussion. Therefore, we will not determine $f$ explicitly.

$F(m,n,p,q)=$$\uparrow/\downarrow$$\uparrow/\downarrow$$E_{n}$$E_{p}$$E_{m}$$E_{q}$$+$$\uparrow/\downarrow$$\uparrow/\downarrow$$E_{n}$$E_{q}$$E_{m}$$E_{p}$$+$$\uparrow/\downarrow$$\uparrow/\downarrow$$E_{n}$$E_{q}$$E_{m}$$E_{p}$
Figure 2: The first two diagrams are totally disconnected. The third is partly
connected and represents the $O(\frac{1}{N^{2}})$ term in $F(m,n,p,q)$.

The function $G$ can be read off from Figure 3. Using the results for single-
trace errors in Eq. (39), we can show that,

$\begin{split}&G_{\uparrow}(m,n,p,q)=(\delta_{m}^{p}\delta_{n}^{q}+\delta_{m}^{q}\delta_{n}^{p})(3n+3m)\\\
&G_{\downarrow}(m,n,p,q)=(\delta_{m}^{p}\delta_{n}^{q}+\delta_{m}^{q}\delta_{n}^{p})(2n+2m).\end{split}$
(42)

From this explicit form, we see that there is indeed a difference between the
these functions, so the logical states are connected by multi-trace operators,
again with $1/N$ suppressed amplitudes provided the indices $m,n,p,q$ are not
too large.

$G_{\uparrow/\downarrow}(m,n,p,q)=$$\uparrow/\downarrow$$E_{p}$$E_{n}$$\uparrow/\downarrow$$E_{m}$$E_{q}$$+$$\uparrow/\downarrow$$E_{q}$$\uparrow/\downarrow$$E_{m}$$E_{n}$$E_{p}$$+$$p\leftrightarrow
q$ Figure 3: The $G$ functions arise from connected diagrams of order
$O(\frac{1}{N^{2}})$ that involve the state
$\operatorname{Tr}(a_{1}^{\dagger})|0\rangle$.

Now, although we only explicitly calculated a particular class of multi-trace
errors, the algorithm generalizes directly to more general cases. The analog
of the $G$ function always consists of connected scattering diagram of logical
states and single-trace operators. In other words, the uncorrectable part of
any multi-trace error can be written as a sum over contributions of the
single-trace errors, due to the factorization property of large $N$ matrix
models (Eq (17)). Therefore, if the aKL condition is satisfied for the single-
trace operators, it will also hold for these multi-trace errors. Note that if
a diagram has four legs connected, then it must be suppressed by
$\frac{1}{N^{2}}$. In special cases we also have three-leg diagrams, where the
error operators involve $\frac{\operatorname{Tr}(a_{1,2}^{\dagger
2})}{\sqrt{2}N}$ or $\frac{\operatorname{Tr}(a_{1,2}^{\dagger})}{\sqrt{N}}$.
For example, the $G$ function in the correlator
$\langle\tilde{\uparrow}|\frac{\operatorname{Tr}(a_{1}^{2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{2}^{2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{1}^{\dagger
2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{2}^{\dagger
2})}{\sqrt{2}N}|\tilde{\uparrow}\rangle$ contains the product of two three-leg
diagrams, as shown in Figure 4. Each of them is proportional to $\frac{1}{N}$,
so the term involving $G$ that depends on the logical state is still
suppressed by $\frac{1}{N^{2}}$.

$\langle\tilde{\uparrow}|\frac{\operatorname{Tr}(a_{1}^{2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{2}^{2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{1}^{\dagger
2})}{\sqrt{2}N}\frac{\operatorname{Tr}(a_{2}^{\dagger
2})}{\sqrt{2}N}|\tilde{\uparrow}\rangle$$\supset$$\uparrow$$\operatorname{Tr}(a_{1}^{2})$$\operatorname{Tr}(a_{2}^{2})$$\uparrow$$\operatorname{Tr}(a_{1}^{\dagger
2})$$\operatorname{Tr}(a_{2}^{\dagger 2})$ Figure 4: Contribution of connected
diagram with three legs. Each is suppressed by $\frac{1}{N}$.

#### Mutual information diagnostic

Now we couple this quantum memory to a single-qubit reference system $R$, and
calculate the mutual information between them (Eq. (29)). At $t=0$, $\rho_{0}$
is a pure state,

$\rho_{0}=\frac{1}{2}\sum_{i,j=1}^{2}|i\tilde{i}\rangle_{SR}\langle
j\tilde{j}|.$ (43)

The system and reference are maximally entangled, so the initial mutual
information is $I(S:R)(0)=2\ln{2}$.

Next we count the effective number of error operators which will contribute
non-trivially to logical errors under the Lindblad dynamics. Operators that
contain more than two $a_{1}$ or $a_{2}$ trivially annihilate the logical
states, so we only need to count the error operators with $r_{2},s_{2}\leq 2$.
The number of gauge invariant error operators with energy $n\omega$ that
satisfy this requirement is bounded by $n^{2}p(n)$, where $p(n)$ is the number
of integer partitions of $n$. So the number of effective errors is sub-
exponential asymptotically.

By construction, the bath model has uniform coupling as required by condition
2 in theorem 1. Together with the aKL condition above, we apply Theorem 1 and
conclude the mutual information decays in the following way,

$I(S:R)(t)=2\ln{2}-K(\frac{t}{N^{2}}),$ (44)

where $K(t)$ is linear in $t$ at early time. More specifically, when $\lambda
t\ll N^{2}$, $K(t)\sim\frac{\sum_{n}n\lambda_{n}t}{N^{2}}$ . It develops a
dependence on higher powers of $t$ at late time $\lambda t\sim N^{2}$, where
$K(t)\sim\sum_{p}(\frac{\sum_{n}n\lambda_{n}t}{N^{2}})^{p}$.

Details of the proof for Eq. (44) and more rigorous results can be found in
Appendix A, but here we briefly sketch the general ideas behind the proof. The
dominant errors in the early-time regime are the single trace errors. In
contrast, in the late time regime, the dominant errors are multi-trace. In
fact, the average size $k$ of errors, defined as the number of traces in an
error operator, grows linearly with time. These errors, when acting on logical
states, cause a loss of quantum information proportional to
$\frac{k^{p}}{N^{2p}}$, leading to the $(\frac{t}{N^{2}})^{p}$ growth behavior
at late time.

Note that the rate of information loss here is captured by the sum over
$\lambda_{n}$, which can be approximated by

$\begin{split}\sum_{n}n^{q}\lambda_{n}=&\sum_{n}n^{q}d_{n}e^{-n\beta\omega}\\\
\sim&\sum_{n}n^{q+2}p(n)e^{-n\beta\omega},\end{split}$ (45)

where $q=1$, and $d_{n}$ is the number of distinct states at the $n$th level.
In this model, $d_{n}\sim n^{2}p(n)$, which grows sub-exponentially with $n$.
As the sum is always finite, it follows that the memory time scales like
$t_{\text{mem}}\sim N^{2}$. From the above discussion, we conclude that this
model serves as a good passive quantum memory.

Now we discuss the zero temperature limit. Since the thermal fluctuation is
suppressed, the relevant error operators are only those which lower the energy
when acting on logical states. In this model, they are
$\frac{\operatorname{Tr}(a_{1,2})}{\sqrt{N}}$ and
$\frac{\operatorname{Tr}(a_{1,2}^{2})}{N}$, whose aKL conditions are given in
Appendix B.1. They lead to a spontaneous decay of quantum information at a
constant rate of $O(\frac{1}{N^{2}})$ at zero temperature. There are also
operators that don’t change energy of the system, but they have no
contribution to the information loss since their thermal spectral functions
(Eq. (22)) vanish at zero temperature limit,
$lim_{\beta\rightarrow\infty}\gamma(\nu=0)=0$.

To modify the encoded information, we need to apply logical operators, which
involve $\operatorname{Tr}(a_{1}^{\dagger}a_{2}^{\dagger})$ or
$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger 2})$. This requires us to
perform non-local operations (or a series of local operations) by bringing the
two modes together and turning on some local coupling, as shown in Figure 5.
We assume this procedure to be accurate. After the logical operations, we can
separate the two modes so that the system again behaves like a passive memory
that protects the encoded information.

$\operatorname{Tr}(a_{1}^{\dagger 2}a_{2}^{\dagger 2})$ (a) Logical operators

$\operatorname{Tr}(a_{1}^{\dagger r})$$\operatorname{Tr}(a_{2}^{\dagger s})$
(b) Thermal errors

Figure 5: In (a), two modes are brought together to implement a logical
operator. In (b), two modes are separated far away. Thermal fluctuation are
sufficiently local so the quantum information is preserved.

We can also implement logical operations without bringing the two modes
together. Instead, we may induce a coupling mediated by a third ancillary mode
$c$. The procedure is that we first bring the mode $c$ close enough to the
mode $a_{1}$ and implement the operator $\operatorname{Tr}(a_{1}^{\dagger
m}c^{\dagger})$. Then we bring together the intermediate mode and the mode
$a_{2}$ and implement the operator $\operatorname{Tr}(a_{2}^{\dagger n}c)$.
These two moves can be combined together to give the desired logical
operation,

$\operatorname{Tr}(a_{2}^{\dagger 2}c)\operatorname{Tr}(a_{1}^{\dagger
2}c^{\dagger})|0\rangle\propto\operatorname{Tr}(a_{1}^{\dagger
2}a_{2}^{\dagger 2})|0\rangle.$ (46)

In summary, this model achieves the power law scaling of memory time with
system size. However, it has some drawbacks. First, it requires a gauge
symmetry that acts consistently on the two matrix modes separated by a large
distance. Second, the logical operation is non-local, and thus can be
difficult to implement. For these reasons, we construct another model, where
the two matrix modes are put together, but with the added price that a high
energy penalty is required.

###  3.2 Local model

In this model, we still consider two matrix harmonic oscillators, but we put
them together and turn on a large coupling between them. The full Hamiltonian
is

$H=\omega\operatorname{Tr}(a_{1}^{\dagger}a_{1})+\omega\operatorname{Tr}(a_{2}^{\dagger}a_{2})+J[\frac{\operatorname{Tr}(a_{1}^{\dagger})\operatorname{Tr}(a_{1})}{N}+\frac{\operatorname{Tr}(a_{2}^{\dagger})\operatorname{Tr}(a_{2})}{N}-1]^{2}+H_{\text{thermal}}$
(47)

with $J\sim 2\omega\ln{N}$. The ground space of this Hamiltonian is spanned by
the following two states,

$\begin{split}|\tilde{\uparrow}\rangle&=\frac{\operatorname{Tr}(a_{1}^{\dagger})}{\sqrt{N}}|0\rangle_{12}\\\
|\tilde{\downarrow}\rangle&=\frac{\operatorname{Tr}(a_{2}^{\dagger})}{\sqrt{N}}|0\rangle_{12}.\end{split}$
(48)

where $|0\rangle_{12}$ is again the tensor product of the vacuum states of the
two oscillators. These states span the code subspace of the model.

We assume that $H_{\rm thermal}$ contains only single trace operators:

$H_{\rm
thermal}=\sum_{\\{n_{k}\\}}\lambda_{\\{n_{k}\\}}b_{\\{n_{k}\\}}\frac{:\operatorname{Tr}(P\\{a_{1}^{\dagger
n_{1}}a_{2}^{\dagger
n_{2}}a_{1}^{n_{3}}a_{2}^{n_{4}}\\}):}{N^{\frac{\sum_{k=1}^{4}n_{k}}{2}}}+h.c..$
(49)

where $P$ again denotes the sum over distinct arrangements of the operators
inside the curly bracket, and all the operators are then normal ordered. In
general, $\operatorname{Tr}(P\\{a_{1}^{\dagger n_{1}}a_{2}^{\dagger
n_{2}}a_{1}^{n_{3}}a_{2}^{n_{4}}\\})$ is a sum over
$\frac{(n_{1}+n_{2}+n_{3}+n_{4})!}{n_{1}!n_{2}!n_{3}!n_{4}!}$ terms.

The interaction term in the Hamiltonian Eq. (47) incurs an energy penalty for
any logical errors that involve $\operatorname{Tr}(a_{1}^{\dagger})$ or
$\operatorname{Tr}(a_{2}^{\dagger})$, so their probability to appear are
suppressed due to the thermal factor. For example, the probability of the
error operator $\operatorname{Tr}(a_{1}^{\dagger})\operatorname{Tr}(a_{2})$ is
$\gamma(J)\gamma(-J)\sim e^{-\frac{J}{T}}$. We set the energy scale $J$ to be
2$\omega\ln{N}$, which guarantees that these uncorrectable errors are unlikely
to occur when $N$ is large enough. All the other operators generated by the
coupling in Eq. (49) are treated as possible errors that could act on the
logical state. To show that these errors are approximately correctable, and to
invoke Theorem 1, it suffices to check that aKL holds for single-trace errors,
$E_{\\{n_{k}\\}}:=\frac{:\operatorname{Tr}(P\\{a_{1}^{\dagger
n_{1}}a_{2}^{\dagger
n_{2}}a_{1}^{n_{3}}a_{2}^{n_{4}}\\}):}{N^{\frac{\sum_{k=1}^{4}n_{k}}{2}}}$,

$\begin{split}\langle\tilde{i}|E_{\\{n_{k}^{\prime}\\}}^{\dagger}E_{\\{n_{k}\\}}|\tilde{j}\rangle&=f_{\\{n_{k}^{\prime}\\},\\{n_{k}\\}}\delta_{ij}+\frac{g_{\\{n_{k}^{\prime}\\},\\{n_{k}\\}}^{ij}}{N^{2}}+O(\frac{1}{N^{4}})\\\
\langle\tilde{i}|E_{\\{n_{k}\\}}|\tilde{j}\rangle&=\frac{e^{ij}_{\\{n_{k}\\}}}{N}+O(\frac{1}{N^{3}}).\end{split}$
(50)

Applying the argument we discussed in Section 3.1, one can again derive the
aKL equation for general multi-trace operators using the above relations for
single-trace operators. This guarantees that the aKL condition holds for all
the error operators generated in this model. A proof of the aKL condition for
the most general error is given in Appendix B.2. Here we discuss a few helpful
examples to show that these errors indeed conform with the aKL condition.

Again, let us first check operators of the form
$\frac{\operatorname{Tr}(a_{1,2}^{\dagger}a_{1,2})}{N}$. They induce a phase
error of size $O(\frac{1}{N})$,

$\begin{split}&\langle\tilde{\uparrow}|\frac{\operatorname{Tr}(a_{1}^{\dagger}a_{1})}{N}|\tilde{\uparrow}\rangle=\frac{1}{N}\\\
&\langle\tilde{\downarrow}|\frac{\operatorname{Tr}(a_{1}^{\dagger}a_{1})}{N}|\tilde{\downarrow}\rangle=0,\end{split}$
(51)

which is consistent with the aKL condition in Eq. (50). Next, we check the aKL
conditions for single-trace errors of the form
$E_{(n,m)}=\frac{\operatorname{Tr}(a_{1}^{\dagger n}a_{2}^{\dagger
m})}{N^{\frac{m+n}{2}}}$, with $n\geq 1,m\geq 1$. For the bit flip errors, we
can show that

$\langle\tilde{\downarrow}|E_{(n^{\prime}+1,m’-1)}^{\dagger}E_{(n,m)}|\tilde{\uparrow}\rangle=\delta_{nn^{\prime}}\delta_{mm^{\prime}}\frac{(n+1)m}{N^{2}}+O(\frac{1}{N^{4}}),$
(52)

and for the phase errors, we have

$\begin{split}&\langle\tilde{\uparrow}|E_{(n^{\prime},m^{\prime})}^{\dagger}E_{(n,m)}|\tilde{\uparrow}\rangle=\delta_{nn^{\prime}}\delta_{mm^{\prime}}[1+\frac{n^{2}-m^{2}}{N^{2}}+O(\frac{1}{N^{4}})]\\\
&\langle\tilde{\downarrow}|E_{(n^{\prime},m^{\prime})}^{\dagger}E_{(n,m)}|\tilde{\downarrow}\rangle=\delta_{nn^{\prime}}\delta_{mm^{\prime}}[1+\frac{m^{2}-n^{2}}{N^{2}}+O(\frac{1}{N^{4}})].\end{split}$
(53)

When $m$ or $n$ equals to $0$, the equations are slightly different. For
instance, when $m=0,n\geq 2$, the matrix elements are,

$\begin{split}&\langle\tilde{\uparrow}|E_{(n^{\prime},0)}^{\dagger}E_{(n,0)}|\tilde{\uparrow}\rangle=\delta_{nn^{\prime}}[1+\frac{n(n-1)}{N^{2}}+O(\frac{1}{N^{4}})]\\\
&\langle\tilde{\downarrow}|E_{(n^{\prime},0)}^{\dagger}E_{(n,0)}|\tilde{\downarrow}\rangle=\delta_{nn^{\prime}}\\\
&\langle\tilde{\uparrow}|E_{(n^{\prime},0)}^{\dagger}E_{(n,0)}|\tilde{\downarrow}\rangle=0.\end{split}$
(54)

Although $H_{\rm thermal}$ only contains single-trace operators, the expansion
in Eq. (28) involves multi-trace operators acting on the state. To see the
effect of these multi-trace errors, we evaluate an example:

$\begin{split}&\langle\tilde{\uparrow}|E_{(p,0)}^{\dagger}E_{(q,0)}^{\dagger}E_{(m,0)}E_{(n,0)}|\tilde{\uparrow}\rangle=F(m,n,p,q)+\frac{G(m,n,p,q)}{N^{2}}+O(\frac{1}{N^{4}})]\\\
&\langle\tilde{\downarrow}|E_{(p,0)}^{\dagger}E_{(q,0)}^{\dagger}E_{(m,0)}E_{(n,0)}|\tilde{\downarrow}\rangle=F(m,n,p,q)+O(\frac{1}{N^{4}}).\end{split}$
(55)

If the Knill-Laflamme conditions are satisfied, then the above expressions
should yield matrix elements that only depend on the errors. We see that the
presence of the $G$ function does allow one to distinguish the logical states,
but the effect is suppressed by $1/N^{2}$. These correlation functions both
have simple diagrammatic representations similar to those shown in Figure 2
and Figure 3, but the explicit form of these two functions are complicated,

$\begin{split}&F(m,n,p,q)=(\delta_{m}^{p}\delta_{n}^{q}+\delta_{m}^{q}\delta_{n}^{p})+\delta_{m+n}^{p+q}\frac{f(p,q,m,n)}{N^{2}}\\\
&G(m,n,p,q)=(\delta_{m}^{p}\delta_{n}^{q}+\delta_{m}^{q}\delta_{n}^{p})[n(n-1)+m(m-1)]+(\delta_{m\pm
1}^{p}\delta_{n\mp 1}^{q}+\delta_{m\pm 1}^{q}\delta_{n\mp
1}^{p})\sqrt{pqmn}.\end{split}$ (56)

with some function $f$ whose exact form is not essential for our argument.

Again all the errors discussed above are approximately correctable errors,
because they can only corrupt the logical information up to
$O(\frac{1}{N^{2}})$. The uncorrectable errors or logical operators are those
that contain $\operatorname{Tr}(a_{1}^{\dagger})$ and
$\operatorname{Tr}(a_{2}^{\dagger})$. However, recall that we have tuned the
coupling $J$ to be proportional to $2\omega\ln{N}$ so that these uncorrectable
errors have a similar $1/N$ suppressed effect comparable in size to the
approximately correctable errors.

By construction, the bath model has uniform coupling as required by condition
2 of Theorem 1. In addition, the number of gauge invariant error operators
with energy $n\omega$ that are induced by Eq. (49) grows asymptotically as
$2^{n}$. In particular, the number of effective error operators that act non-
trivially on the logical state, i.e., ones contributing to the information
leak in the Lindblad equation, must have $n_{2},n_{4}\leq 1$. Therefore, the
total number of effective errors with energy $\omega n$ is bounded by an
exponential function, as required by condition 1.

Because the aKL condition holds and the bath model only generates errors which
are protected by the code, we again invoke the Theorem 1. The mutual
information between the encoded qubit $S$ and a reference $R$ is

$I(S:R)(t)=2\ln{2}-K_{\rm local}(\frac{t}{N^{2}}).$ (57)

The expression $K_{\rm local}(t)$ in this model is slightly different from
$K(t)$, since the $O(\frac{1}{N^{2}})$ term in the aKL condition depends
quadratically on $n$ (see Eq. (52)). When $\lambda t\ll N^{2}$, $K_{\rm
local}(\frac{t}{N^{2}})\sim\frac{\sum_{n}n^{2}\lambda_{n}t}{N^{2}}$. At late
time $\lambda t\sim N^{2}$, $K_{\rm
local}(t)\sim(\sum_{n}n^{2}\lambda_{n}t)^{2}$.

More importantly, note that the rate of information loss is given by

$\begin{split}\sum_{n}n^{q}\lambda_{n}=&\sum_{n}n^{q}d_{n}e^{-n\beta\omega}\\\
\sim&\sum_{n}n^{q}2^{n}e^{-n\beta\omega},\end{split}$ (58)

where $q=2$, and $d_{n}\sim 2^{n}$ is the number of single trace states at the
$n$th energy level. A crucial distinction from the previous model is the
Hagedorn-like transition with a critical temperature
$T_{c}=\frac{\omega}{\ln{2}}$ [26]. In the low temperature phase where
$T<T_{c}$, the rate of information loss (58) is finite and $N$-independent.
Therefore the large $N$ gauged matrix model has a long memory time
$t_{\text{mem}}\sim N^{2}$. However, for $T>T_{c}$, (58) exhibits a Hagedorn-
like behaviour888Note that for finite $N$, the Hagedorn behaviour does not
continue to arbitrary energies, but only to order $N^{2}$, the total number of
gauge invariant degrees of freedom. The counting of gauge invariant operators
can be found in Section 2.4. . As a result, $K_{\rm local}(t)$ now also scales
with $N$. This $N$ dependence cancels out the large $N$ suppression, and the
memory time no longer scales with the system size. This behaviour is similar
to the holographic code [10], where the temperature needs to be lower than a
critical value to maintain the robustness of the encoded information.

Since the logical states are within the ground space of this model, there are
no operators that can lower the energy. The zero energy error operators do not
contribute to information loss at zero temperature limit, due to the vanishing
thermal spectral function $\lim_{\beta\rightarrow\infty}\gamma(\nu=0)$.
Therefore, when $\beta\rightarrow\infty$, we expect the information decay rate
to scale as $\frac{1}{N^{2}}e^{-\beta\omega}$ according to Eq. (58).

##  4 Relaxing the gauge constraint

We now relax the gauge constraint to a global symmetry, which we then impose
energetically. The purpose of this is two-fold — (i) we want to understand how
the stability of quantum memory is tied to the gauge symmetry specifically,
and (ii) there are very few models in nature that are known to have gauge
symmetries. For example, Quantum ElectroDynamics (QED) has a $U(1)$ symmetry
and Quantum ChromoDynamics (QCD) has an $SU(3)$ symmetry. However, there is no
known $SU(N)$ gauged model that can be realized in the laboratory, especially
for arbitrarily large $N$. Therefore, in this section, we consider a single
oscillator model with only a global $SU(N)$ symmetry. The system of interest
has Hamiltonian

$\begin{split}H_{\rm
Gauss}=\frac{\mathcal{J}}{N}\sum_{ij}G^{i}_{j}G^{j}_{i}\\\
G^{i}_{j}=\sum_{k}a^{\dagger i}_{k}a^{k}_{j}-a^{\dagger
k}_{j}a^{i}_{k},\end{split}$ (59)

where $G^{i}_{j}|\psi\rangle=0$ is the Gauss law constraint associated with
the gauge symmetry. The energy penalty term plays the role of a Casimir
operator in the group $SU(N)$, with eigenvalues proportional to $\mathcal{J}$
and being the same for all states in an irreducible representation. This term
in the Hamiltonian forces the low energy states to be singlets of the symmetry
group.

In the absence of gauge symmetry, we allow the modes $a^{i}_{j}$’s to couple
independently with the thermal bath operator. Let us consider the following
thermal coupling Hamiltonian,

$\begin{split}H_{\rm
thermal}=&\lambda_{1}\sum_{ij}b^{i}_{j}a^{j(\dagger)}_{i}+\lambda_{2}\sum_{ijkl}b^{ik}_{jl}a^{j(\dagger)}_{i}a^{l(\dagger)}_{k}+\lambda_{3}\sum_{ijklpq}b^{ikp}_{jlq}a^{j(\dagger)}_{i}a^{l(\dagger)}_{k}a^{q(\dagger)}_{p}+\cdots\\\
&+h.c..\end{split}$ (60)

One can show that the Lindblad Master equation takes the form of

$\dot{\rho}(t)=\lambda_{1}\frac{\operatorname{Tr}(a^{\dagger})}{\sqrt{N}}\rho\frac{\operatorname{Tr}(a)}{\sqrt{N}}+\lambda_{1}e^{-\frac{\mathcal{J}}{T}}\sum_{ij}\tilde{a}^{\dagger
i}_{j}\rho\tilde{a}^{j}_{i}+\cdots,$ (61)

where $\tilde{a}^{i}_{j}$ represents the component that transforms in the
adjoint representation of $SU(N)$, such that

$\tilde{a}^{\dagger i}_{j}=a^{\dagger
i}_{j}-\frac{\operatorname{Tr}(a^{\dagger})}{N}\delta^{i}_{j}.$ (62)

We show in Appendix C that it is a consequence of the $H_{\rm Gauss}$ term in
the Hamiltonian.

These non-singlet errors acting on the “gauge invariant” logical states, i.e.,
logical states satisfying the Gauss law constraint, only corrupt quantum
information at $O(\frac{1}{N^{2}})$, which we show in Appendix B.3. The
summation in Eq. (61) is over $N^{2}$ number of such errors and therefore the
total effect of these non-singlet errors is proportional to
$N^{2}e^{-\frac{\mathcal{J}}{T}}$. Hence the coupling strength $\mathcal{J}$
has to be at least $\sim 2T\log{N}$ in order to suppress their effects to the
same order as those from the the singlet errors. Note that this model
effectively contains one such matrix oscillator and does not exhibit a
Hagedorn-like transition. Nevertheless, it is clear that stability is only
ensured when $T\lesssim\mathcal{J}$, hence can also depend on the system
size999This should be contrasted with the models with exact gauge symmetry,
where the temperature that sets apart a stable and an unstable memory only
depends on the natural frequency of the system. Note that in the local model
with gauge symmetry, there is an additional scale introduced by $J$ which has
to satisfy $J\gtrsim\omega\log N$ to to ensure stability. However, as long as
we have tuned $J>T_{c}$, then $T_{c}$ will be the relevant energy scale for
the low energy stability.. In the following, we will assume their suppression
with a strong enough coupling. As for the singlet error operators, they are
defined as,

$E_{(\\{n_{k},m_{k})\\}}:=\prod_{k}\frac{\operatorname{Tr}(P\\{a^{\dagger
n_{k}}a^{m_{k}}\\})}{N^{\frac{\sum_{k}n_{k}+m_{k}}{2}}}.$ (63)

The logical states we choose in this model are (divided by proper
normalization constant),

$\begin{split}&|\tilde{1}\rangle=\frac{\operatorname{Tr}(a^{\dagger
L})^{2}}{\mathcal{N}(L)}|0\rangle\\\
&|\tilde{2}\rangle=\frac{\operatorname{Tr}(a^{\dagger
L+1})\operatorname{Tr}(a^{\dagger L-1})}{\mathcal{N}(L)}|0\rangle.\end{split}$
(64)

We also define the orthogonormal basis,

$\begin{split}&|\tilde{\uparrow}\rangle=\frac{|\tilde{1}\rangle+|\tilde{2}\rangle}{2(1+\langle\tilde{1}|\tilde{2}\rangle)}\\\
&|\tilde{\downarrow}\rangle=\frac{|\tilde{1}\rangle-|\tilde{2}\rangle}{2(1-\langle\tilde{1}|\tilde{2}\rangle)}.\end{split}$
(65)

In order to recast the code subspace as a lower energy subspace, let us
consider a shifted version of the system Hamiltonian and write the full
Hamiltonian as

$H=\frac{\mathcal{J}}{N}\sum_{ij}G^{i}_{j}G^{j}_{i}+\omega(\operatorname{Tr}(a^{\dagger}a)-L)^{2}+H_{\rm
thermal},$ (66)

where $L>0$ is some arbitrary integer.

One can check that the aKL condition,

$\begin{split}\langle\tilde{i}|E_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\}}^{\dagger}E_{\\{(n_{k},m_{k})\\}}|\tilde{j}\rangle&=f_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\},\\{(n_{k},m_{k})\\}}\delta_{ij}+\frac{g_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\},\\{(n_{k},m_{k})\\}}^{ij}}{N^{2}}+O(\frac{1}{N^{4}})\\\
\langle\tilde{i}|E_{\\{(n_{k},m_{k})\\}}|\tilde{j}\rangle&=\frac{e^{ij}_{\\{(n_{k},m_{k})\\}}}{N}+O(\frac{1}{N^{3}}),\end{split}$
(67)

holds for all relevant errors induced by (63) as long as they do not contain
$\operatorname{Tr}(a^{\dagger L})$ or $\operatorname{Tr}(a^{\dagger L\pm 1})$
as a factor. For example, consider the class of error operators
$E_{n}=\frac{\operatorname{Tr}(a^{\dagger n})}{\sqrt{n}N^{\frac{n}{2}}}$. The
relevant matrix elements are

$\begin{split}&\langle\tilde{\uparrow}|E_{m}^{\dagger}E_{n}|\tilde{\downarrow}\rangle=\delta_{mn}[\frac{2n}{N^{2}}+O(\frac{1}{N^{4}})]\\\
&\langle\tilde{\uparrow}|E_{m}^{\dagger}E_{n}|\tilde{\uparrow}\rangle-\langle\tilde{\downarrow}|E_{m}^{\dagger}E_{n}|\tilde{\downarrow}\rangle=O(\frac{1}{N^{4}}).\end{split}$
(68)

Note that although the Hamiltonian and logical states depend on $L$, the error
induced by a low energy operator has no such dependence (the
$O(\frac{1}{N^{2}})$ term does not depend on $L$); ergo, the aKL conditon in
the standard form Eq. (18) holds for such singlet errors (63). We present the
details of the proof in Appendix B.3).

Similar to the discussion in Section 3.1, the number of singlet errors with
the energy $n^{2}\omega$ scales as the number of integer partitions $p(n)$.
This is clearly slower than exponential growth. Therefore, if the thermal
coupling in (60) does not contain $\operatorname{Tr}(a^{\dagger L})$ and
$\operatorname{Tr}(a^{\dagger L\pm 1})$ as factors, and $\lambda\sim O(1)$,
then Theorem 1 applies and we again attain a memory time quadratic in $N$.
Since the logical states are within the ground space, there are no spontaneous
decay of quantum information. So the memory time is proportional to
$N^{2}e^{\frac{\omega}{T}}$ at zero temperature limit.

However, the low energy errors like $\operatorname{Tr}(a^{\dagger
L})\operatorname{Tr}(a^{L})$ do appear, which means that the aKL condition
does not hold for these operators. Nevertheless, it is possible to retain this
long memory time in this model as long as we also suppress the singlet errors
containing such terms. One way to achieve this is to assume that the coupling
constants are exponentially suppressed with respect to the number of modes
involved, i.e., $\lambda_{k}\sim\lambda^{k}$ for some $\lambda<1$. This gaps
out the dangerous error for large $L\sim\log{N}$ by suppressing the error
probability to order $\frac{1}{N^{2}}$. In the meantime, it is also a natural
requirement in this model because the thermal bath mode couples to each
oscillator individually.

One potential drawback in this model (and also the fully gauged model in
Section 3.2) is that the coupling constants $\mathcal{J},J\propto\ln N$ depend
on the system size. However, this is not as unreasonable as it looks because
the interaction energy per harmonic oscillator mode is small101010The energy
for thermal state at temperature $T\sim\omega\ll\mathcal{J}$ is $O(1)$ times
$\mathcal{J}$. See Appendix D for more detail.  and is proportional to
$\frac{\ln{N}}{N^{2}}\omega$. For reference, let us contrast this
$N$-dependence with that of the toric code model with system size equal to
$N$. Recall that the toric code is a $Z_{2}$ lattice gauge theory and the
Hamiltonian also involves terms that are proportional to the Gauss’s law
operators, with coupling constant $\Delta$. When it is coupled to a bath with
fixed temperature $T$, the memory time is proportional to
$\frac{1}{N}e^{\frac{\Delta}{T}}$. It becomes a viable quantum memory if one
can tune the coupling $\Delta$ to $T\ln{N}$, which corresponds to an energy
density much larger than that required in our matrix model.

##  5 Large N spin model

In the previous sections, we saw that exact $SU(N)$ symmetries in oscillator
models can lead to robust quantum memories, both when the symmetry is fully
gauged and when the singlet restriction is enforced energetically. It is
natural to ask whether such exact symmetries and infinite dimensional Hilbert
spaces are truly necessary for robustness. In other words, can we realize an
analogous phenomenon with finite dimensional systems, e.g. involving qubits,
that only have an approximate $SU(N)$ symmetry? In this section, we propose
one such model and show that the requisite large $N$ factorization properties
are mostly preserved up to small corrections in the low energy sector. As a
result, the general conclusions of Theorem 1 still apply to this system in a
way that is similar to Section 4.

To motivate the model, let us first review how one can approximately embed the
low-energy states of an oscillator into a spin Hilbert space. Consider a
spin-$j$ representation of $SU(2)$. When acting on states with $S^{z}$
eigenvalue close to maximal, $m\lesssim j$, the spin commutator becomes

$[S^{x},S^{y}]=i\hbar S^{z}\approx i\hbar^{2}j+\cdots.$ (69)

This suggest that $S^{x}$ and $S^{y}$ can be approximately mapped to a
canonical pair. Setting

$Q=\frac{\ell}{\hbar\sqrt{j}}S^{x},\,\,\,\,P=\frac{1}{\ell\sqrt{j}}S^{y}$ (70)

for some as yet undetermined length $\ell$, we see that

$[Q,P]\approx i\hbar+\cdots.$ (71)

Now considering the harmonic oscillator Hamiltonian

$H=\frac{1}{2}(P^{2}+\omega^{2}Q^{2})=\frac{1}{2}\left(\frac{1}{\ell^{2}j}(S^{y})^{2}+\frac{\omega^{2}\ell^{2}}{\hbar^{2}j}(S^{x})^{2}\right),$
(72)

we may set $\ell^{2}=\hbar/\omega$ to obtain

$H=\frac{\omega}{2\hbar
j}\left((S^{x})^{2}+(S^{y})^{2}\right)=\frac{\omega}{2\hbar
j}(S^{2}-(S^{z})^{2}).$ (73)

Using the standard spectrum of $S^{2}$ and $S^{z}$ with $m=j-\Delta$, the
spectrum of $H$ is

$\frac{\hbar\omega}{2j}\left(j(j+1)-(j-\Delta)^{2}\right)=\hbar\omega\left(\Delta+\frac{1}{2}+O(\Delta^{2}/j)\right).$
(74)

From this discussion, we expect that the physics of our oscillator models may
be suitably extended by replacing the oscillators with qudits of dimension $d$
($d=2j+1$ in the above mapping). At large $d$ we have a nearly identical low-
energy spectrum, but the models make sense for any $d$ and we will show that
even $d=2$ ($j=1/2$) suffices.

Specifically, we can construct a system of $N^{2}$ spins, which can be thought
of as qudits with local dimension $d$. It is tempting to write each spin
operator in the adjoint representation of $SU(N)$, meaning each spin operator
carries indices $i,j$ running from $1$ to $N$ and transforms as

$S^{i\pm}_{j}\rightarrow\sum_{kl}U^{i}_{k}S^{k\pm}_{l}U^{\dagger l}_{j}.$ (75)

Each spin also satisfies the following $SU(2)$ algebra,

$\begin{split}&[S_{i}^{j+},S_{j}^{i-}]=2S_{i}^{jz}\\\
&[S_{i}^{jz},S_{i}^{j+}]=S_{i}^{j+}\\\
&[S_{i}^{jz},S_{j}^{i-}]=-S_{j}^{i-},\end{split}$ (76)

where $i,j$ again label the different spins or qudits. However, the $SU(N)$
symmetry does not act consistently on the spin algebra due to the first line
above: we cannot interpret $S_{i}^{jz}$ as being in the adjoint of $SU(N)$
while also obeying the commutation relation. Moreover, even if we grant that
$S_{i}^{jz}$ is not in the adjoint representation, it is impossible to
construct $SU(N)$ generators $G_{i}^{j}$ with the desired commutation
relations with $S_{i}^{j\pm}$, i.e. there is no $G_{i}^{j}$ obeying111111The
inconsistency of $SU(N)$ and the $SU(2)$ algebra can be seen from the equation
$\sum_{k}[[S^{i+}_{k},S^{k-}_{j}],S^{l+}_{h}]=\delta^{i}_{j}\delta^{j}_{h}S^{l+}_{h}$,
where the right hand side transforms differently from the left hand side.

$[G_{i}^{j},S_{k}^{l\pm}]=\delta_{i}^{l}S_{k}^{j\pm}-\delta_{k}^{j}S_{i}^{l\pm}.$
(77)

However, it is still possible to realize an approximate symmetry in terms of
the pseudo-generators,

$\begin{split}\tilde{G}_{i}^{j}=\sum_{k}S^{j+}_{k}S^{k-}_{i}-S^{k+}_{i}S^{i-}_{k}\\\
[\tilde{G}_{i}^{j},S_{k}^{l\pm}]=(\delta_{i}^{l}S_{k}^{j\pm}-\delta_{k}^{j}S_{i}^{l\pm})S^{lz}_{k}.\end{split}$
(78)

Indeed, if $S^{z}$ could be approximated as a constant operator, then we would
have the desired commutators. This approximate algebra thus connects back with
the oscillator models of the previous sections where the symmetry becomes
exact. In what follows, we work with the pseudo-generators $\tilde{G}_{i}^{j}$
and set $d=2$.

Similar to our recipe in the previous sections, one can construct the analogs
of the gauge invariant states by acting multi-trace operators on the tensor
product of spin down states

$\begin{split}&\operatorname{Tr}(S^{+n_{1}})\operatorname{Tr}(S^{+n_{2}})\cdots\operatorname{Tr}(S^{+n_{k}})|0\rangle\\\
&|0\rangle:=|\downarrow\downarrow\cdots\downarrow\rangle.\end{split}$ (79)

For convenience, we still refer to these analogues of the gauge invariant
operators (states) singlet operators (states). Indeed, one can see that the
desired symmetry is only approximate because the algebra of such operators is
not closed. For example,

$\operatorname{Tr}(S^{-2})\operatorname{Tr}(S^{+2})\operatorname{Tr}(S^{+2})|0\rangle=(4N^{2}-4N-8)\operatorname{Tr}(S^{+2})|0\rangle+32\sum_{i}S^{i+}_{i}S^{i+}_{i}|0\rangle.$
(80)

Generally, the product of such singlet terms can produce non-singlet terms
that are only sub-leading in $N$ 121212If we further restrict ourselves to
only consider spin-$\frac{1}{2}$,i.e., $d=2$, then the non-singlet terms also
vanishes. This can be a convenient simplification when the local systems are
exactly qubits, however, our analysis of robustness will also hold for systems
with $d>2$..

Now let us consider a system Hamiltonian $H_{S}=hH_{0}+\mathcal{J}H_{G}$,
where

$H_{0}=(\sum_{ij}S^{iz}_{j}-S_{\rm tot})^{2},$ (81)

and $-\frac{N^{2}}{2}\leq S_{\rm tot}\leq\frac{N^{2}}{2}$ is some
$N$-dependent offset. $H_{G}$ is an energy penalty term similar to the Gauss
law constraint; it corresponds to the Casimir operator that enforces the
approximate $SU(N)$ symmetry for low energy states131313The Casimir operator
$C=\sum_{ij}G^{i}_{j}G^{j}_{i}$ for an exact $SU(N)$ symmetry has a spectrum
of the form $NwJ$, where $w$ is the highest weight of given representation.
Analogously, the Hamiltonian $H_{G}$ also has an energy spectrum of the form
of $wJ$ up to $\frac{1}{N}$ corrections. ,

$H_{G}=\frac{1}{N}\tilde{G}^{2}:=\frac{1}{N}\sum_{ij}\tilde{G}_{i}^{j}\tilde{G}_{j}^{i}.$
(82)

When $\mathcal{J}$ is large, it creates a large gap in the spectrum. We define
the low energy sector as the states with eigenvalues suppressed by
$\frac{1}{N}$. The energy eigenstates of $H_{S}$ are generally superpositions
of the singlet and non-singlet states, but it can be shown that the low energy
sector is isomorphic to the singlet subspace. For example, the following are
four states created by product of three $S^{+}$’s, which are simultaneous
eigenstates of both $H_{0}$ and $H_{G}$,

$\begin{split}&|\psi_{1}\rangle=\operatorname{Tr}(S^{+3})|0\rangle\\\
&|\psi_{2}\rangle=[\operatorname{Tr}(S^{+})^{3}-3\operatorname{Tr}(S^{+2})\operatorname{Tr}(S^{+})]|0\rangle\\\
&|\psi_{3}\rangle=[\operatorname{Tr}(S^{+2})\operatorname{Tr}(S^{+})+(1+\frac{2}{N})\operatorname{Tr}(S^{+})^{3}-\frac{19}{2N}\sum_{i}(S^{+2})^{i}_{i}S^{+i}_{i}]|0\rangle\\\
&|\psi_{4}\rangle=[\sum_{i}(S^{+2})^{i}_{i}S^{+i}_{i}-\frac{1}{N}\operatorname{Tr}(S^{+2})\operatorname{Tr}(S^{+})]|0\rangle.\end{split}$
(83)

The first three states are in the low energy sector, since they are
eigenstates of $H_{G}$ with eigenvalues $E_{1}^{G}=E_{2}^{G}=0$ and
$E_{3}^{G}=\frac{7\mathcal{J}}{2N}$. The fourth state has eigenvalue
$E_{4}^{G}=\mathcal{J}$, and hence a high energy excitation.

There are two potential causes for this model to deviate from the one in
Section 4 — effects from the symmetry being approximate and finite
dimensionality of each spin. Formally, we should still expect the low energy
states, i.e., ones with energy $\ll\mathcal{J}$, to have a one-to-one
correspondence with the multi-trace operators formed by $S^{+}$’s. More
precisely, their low energy spectra should coincide up to corrections that are
$1/N$ suppressed and the multiplicity, when we ignore such corrections, should
be identical.

We see this by noting that their energy only depend on number of $S^{+}$
involved, with small corrections proportional to $\frac{\mathcal{J}}{N}$
(Figure 6). Therefore, up to these small corrections introduced by the
deviations from the exact symmetry, the operator counting at each energy level
is the same as that in Section 4. These corrections can be connected to the
break down of the exact symmetry.

As for the correction from finite dimensionality, thanks to the “Gauss law”
constraint $H_{G}$ and the form of $H_{0}$, these corrections would only
significantly modify this counting when the total number of $S^{+}$ involved
is of $O(N)$. Therefore, in the low energy regime where we operate, both
effects are suppressed by having a large enough $N$.

$n=0$$n=1$$n=2$$n=3$$\frac{\mathcal{J}}{N}$ Figure 6: Specturm of the low
energy states. The ground state has energy zero. Excited states has energy
proportional to $n^{2}h$. Note that due to the non-exactness of symmetry,
there are small energy splitting proportional to $\frac{\mathcal{J}}{N}$ . For
sufficiently large $N$, this small gap can be neglected and the degeneracy at
each level is the same as in harmonic oscillator model. Non-singlet
excitations have energy $O(\mathcal{J})$, which we do not count as low energy.

Similar to the harmonic oscillator model, we consider the full Hamiltonian,

$H=hH_{0}+\mathcal{J}H_{G}+H_{\rm thermal}$ (84)

with $H_{\rm thermal}$ being

$\begin{split}H_{\rm
thermal}=&\lambda_{1}\sum_{ij}b^{i}_{j}S^{j\pm}_{i}+\lambda_{2}\sum_{ijkl}b^{ik}_{jl}S^{j\pm}_{i}S^{l\pm}_{k}+\lambda_{3}\sum_{ijklpq}b^{ikp}_{jlq}S^{j\pm}_{i}S^{l\pm}_{k}S^{q\pm}_{p}+\cdots\\\
&+h.c..\end{split}$ (85)

The lowest energy states must have total spin equal to $S_{\rm tot}$. Let us
define $L=\frac{S_{\rm tot}+N^{2}/2}{2}$ and pick a two dimensional subspace
as the code subspace. With proper normalization,

$\begin{split}|\tilde{1}\rangle&=\frac{\operatorname{Tr}(S^{+L})^{2}}{\mathcal{N}}|0\rangle\\\
|\tilde{2}\rangle&=\frac{\operatorname{Tr}(S^{+L+1})\operatorname{Tr}(S^{+L-1})}{\mathcal{N}}|0\rangle.\\\
\end{split}$ (86)

Again, we define the orthogonormal basis in this logical subspace as,

$\begin{split}&|\tilde{\uparrow}\rangle=\frac{|\tilde{1}\rangle+|\tilde{2}\rangle}{2(1+\langle\tilde{1}|\tilde{2}\rangle)}\\\
&|\tilde{\downarrow}\rangle=\frac{|\tilde{1}\rangle-|\tilde{2}\rangle}{2(1-\langle\tilde{1}|\tilde{2}\rangle)}.\end{split}$
(87)

with the overlap
$\langle\tilde{1}|\tilde{2}\rangle=\frac{f(L)L\sqrt{L^{2}-1}}{\sqrt{2}N^{2}}$,
where $f(L)$ is a complicated function defined in Appendix B.4.

Let us couple this system to a thermal bath at temperature $T$. First note
that the non-singlet excitations are again suppressed by a factor of
$e^{-\frac{\mathcal{J}}{T}}$ because of the energy penalty term $H_{G}$, while
the singlet interaction has no such suppression. Then let us reorganize the
system-bath interactions by the approximate symmetry such that

$\begin{split}H_{\rm
thermal}&=\lambda_{1}\sum_{ij}b_{j}^{i}S^{j\pm}_{i}+\lambda_{2}\sum_{ijkl}b^{ik}_{jl}S^{j\pm}_{i}S^{l\pm}_{k}+\cdots+h.c.\\\
&=\lambda_{1}\sum_{ij}b_{j}^{i}(S_{i}^{j+}-\frac{\operatorname{Tr}(S^{+})}{N}\delta_{i}^{j})+\lambda_{1}\frac{\operatorname{Tr}(b)}{\sqrt{N}}\frac{\operatorname{Tr}(S^{+})}{\sqrt{N}}+\cdots.\end{split}$
(88)

Similar to the discussion in last section, this interaction induces an order
$N^{2}$ number of non-singlet errors. Therefore, the total error contribution
from non-singlet excitations is $\sim N^{2}e^{-\frac{\mathcal{J}}{T}}$ times
that of singlet excitations. To suppress both types of errors to the same
order, we require the coupling $\mathcal{J}$ to be greater than $2h\log{N}$ .

More precisely, the singlet error operators are defined as

$E_{(\\{n_{k},m_{k})\\}}:=\prod_{k}\frac{\operatorname{Tr}(P\\{S^{+n_{k}}S^{-m_{k}}\\})}{N^{\frac{\sum_{k}n_{k}+m_{k}}{2}}}.$
(89)

Again, we show that the aKL condition

$\begin{split}\langle\tilde{i}|E_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\}}^{\dagger}E_{\\{(n_{k},m_{k})\\}}|\tilde{j}\rangle&=f_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\},\\{(n_{k},m_{k})\\}}\delta_{ij}+\frac{g_{\\{(n_{k}^{\prime},m_{k}^{\prime})\\},\\{(n_{k},m_{k})\\}}^{ij}}{N^{2}}+O(\frac{1}{N^{4}})\\\
\langle\tilde{i}|E_{\\{(n_{k},m_{k})\\}}|\tilde{j}\rangle&=\frac{e^{ij}_{\\{(n_{k},m_{k})\\}}}{N}+O(\frac{1}{N^{3}})\end{split}$
(90)

is satisfied in Appendix B.4. For example, one can see that the operator
$\frac{\operatorname{Tr}(S^{+}S^{-})}{N}=\frac{\sum_{ij}S^{iz}_{j}}{N}$ acts
as the logical identity, since the logical states have the same total spin in
$z$ direction. For another example, the aKL conditions for errors
$E_{n}=\frac{\operatorname{Tr}(S^{+n})}{\sqrt{n}N^{\frac{n}{2}}}$ are

$\begin{split}&\langle\tilde{\uparrow}|E_{m}^{\dagger}E_{n}|\tilde{\downarrow}\rangle=\delta_{mn}[\frac{2n}{N^{2}}+O(\frac{1}{N^{4}})]\\\
&\langle\tilde{\uparrow}|E_{m}^{\dagger}E_{n}|\tilde{\uparrow}\rangle-\langle\tilde{\downarrow}|E_{m}^{\dagger}E_{n}|\tilde{\downarrow}\rangle=O(\frac{1}{N^{4}}).\end{split}$
(91)

Note that in these examples, and more generally in Appendix B.4, deviations
from the KL condition only depend on $n$ but not on $L$. Therefore, the aKL
condition remains valid even for very large $L$. Similar to the last section,
we can pick $L\gtrsim\log{N}$ to make sure uncorrectable errors are suppressed
for large $N$. For instance, one may consider choosing $L=N^{2}/4$ such that
$S_{\rm tot}=0$.

Since this spin model has a structure very similar to the ungauged harmonic
oscillator model, we can follow the same argument in Section 4. In the low
energy sector with $\epsilon\ll\mathcal{J}\sim\log N$, (approximate) singlet
errors are the dominant error modes. The number of such error operators grows
sub-exponentially with energy because these operators still have the multi-
trace structure whose counting are the same as the one in Section 4. There are
small corrections because of the approximate symmetry, however they are
proportional to $\frac{\mathcal{J}}{N}$, and therefore for large $N$, they do
not affect the overall structure of the spectrum for the purpose of our
analysis (Figure 6). Similarly, we also need the thermal coupling strength to
decay exponentially with number of spins to suppress the singlet error
operators such as $\operatorname{Tr}(S^{-L})\operatorname{Tr}(S^{+L})$. As
argued in the previous section, this is a natural assumption since the bath
mode couples to each spin individually. Since the rest of the argument is
identical to those in Section 4, we can then invoke Theorem 1 to conclude that
the memory time $t_{m}$ scales quadratically with $N$. There is no operator
that can lower the logical states’ energy, so the zero temperature limit of
memory time is still proportional to $N^{2}e^{\frac{\omega}{T}}$. We remark
that this model also has the same drawback as the ungauged harmonic oscillator
model, that is, the interactions remain non-local and we need the coupling
strength $\mathcal{J}$ to scale as at least $\log{N}$.

##  6 Discussion and Conclusion

Motivated by connections between error correction and phenomena in gauge
theory and quantum gravity, we have taken a first step in understanding how
$SU(N)$ gauge symmetries at large $N$ can be used to construct stable quantum
memories. In particular, we construct a gauged matrix quantum mechanics toy
model with non-local logical operations that is a self-correcting quantum
memory. The memory time of this model is polynomial in $N$ when coupled to a
thermal bath at non-zero temperature. The presence of gauge symmetry induces
approximate quantum error correcting codes that become exact in the large $N$
limit, thanks to the large $N$ factorization property. Furthermore, gauge
invariance limits the number of error patterns that can be excited up to a
certain energy scale. The combination of these factors suppresses the
proliferation of memory errors under thermal noise. A similar model with local
logical operation is also possible, but as a trade off, we impose an
additional energy penalty with a coupling of strength $\sim\log N$. The memory
time is again polynomial in $N$, but only for temperatures below a critical
value.

As large $N$ gauge symmetries are not known to occur in nature, we relax the
gauge symmetries and impose only large $N$ global symmetries energetically on
a single matrix oscillator. We also construct a non-local finite dimensional
spin model which further relaxes the symmetry constraint where we only enforce
approximate large $N$ symmetries energetically. Interestingly, these
constructions remain robust in both scenarios under the assumed system-bath
coupling as long as the symmetry constraints can be energetically enforced
with a penalty that scales as $O(\log N)$.

Stepping back from these particular examples, our work provides one set of
guidelines for building large $N$ passive quantum memories in which the
encoded information remains robust when coupled to a thermal bath. Sufficient
conditions on the quantum system and its coupling to the bath are captured in
Theorem 1. We hope to extend our current work in a number of future
directions. Let us roughly categorize these potential generalizations in the
direction of optimality, dimensionality, symmetries, variety, and
practicality.

Although we can construct models that have $\text{poly}(N)$ memory time, they
also often come with other features that are practically difficult, such as
non-local interactions and $\log N$ coupling strengths. However, we by no
means claim that these are the optimal constructions. Indeed, it is possible
that a more careful design can produce better self-correcting memories. For
example, an interesting theoretical problem is whether there exists a gauged
$SU(N)$ oscillator model with local coupling for which the quantum memory is
stable against thermal noise.

It is also natural to generalize these findings to systems with gauge
symmetries in higher dimensions by considering geometrically coupled
oscillators, lattice gauge theories, and concatenations with other quantum
error correcting codes. We will take a first step towards tackling a subset of
these problems in an upcoming work [12]. Although we enforce global symmetries
energetically in our work, one can similarly ask whether there are advantages
in considering systems with emergent gauge symmetries at low energies, which
is better understood in higher dimensional systems. Relatedly, one can examine
what aspects of the $SU(N)$ symmetry are actually essential. For instance,
does a similar result hold for discrete, instead of continuous, non-abelian
symmetries?

Drawing intuition from gauge theories, we can also examine how to generalize
such models by coupling the “pure gauge” degrees of freedom to fundamental
degrees of freedom like matter fields and defects. In particular, how does the
configuration of these other degrees of freedom impact our conclusions? One
may, for instance, introduce a variety of fields that may interact under these
gauge couplings. Extensions in this direction draw heavy parallels with
traditional theoretical particle physics, and can open up new connections with
quantum information in this area.

Lastly, as a practical aspect of any quantum memory, information needs to be
encoded and extracted relatively easily. This is yet unexplored for both the
infinite dimensional oscillator and the finite dimensional spin system. As a
standalone quantum code, the gauged matrix oscillator is a continuous variable
system that resembles bosonic codes [32] in many ways. However, we have barely
begun to analyze its code properties. It is desirable to understand its
similarities and differences with existing bosonic codes, its robustness
against different kinds of errors, and how active error correction like
syndrome extraction can be performed.

## Acknowledgement

We thank Alexey Milehkin, John Preskill for helpful comments and discussions.
C.C. acknowledges the support by the U.S. Department of Defense and NIST
through the Hartree Postdoctoral Fellowship at QuICS, the Air Force Office of
Scientific Research (FA9550-19-1-0360), and the National Science Foundation
(PHY-1733907). G.C. acknowledges support from the U.S. Department of Energy,
Office of Science, Office of Advanced Scientific Computing Research,
Accelerated Research for Quantum Computing program “FAR-QC”. B.G.S.
acknowledges support from the AFOSR under grant number FA9550-19-1-0360.

## References

  * [1] A.Y. Kitaev.  Proceedings of the 3rd International Conference of Quantum Communication and Measurement, Ed. O. Hirota, A. S. Holevo, and C. M. Caves.  Plenum, New York, 1997.
  * [2] A.Yu. Kitaev.  Fault-tolerant quantum computation by anyons.  Annals of Physics, 303(1):2–30, 2003.
  * [3] Alexei Kitaev and John Preskill.  Topological Entanglement Entropy.  Phys. Rev. Lett., 96(11):110404, March 2006.
  * [4] R. Alicki, M. Fannes, and M. Horodecki.  On thermalization in Kitaev’s 2D model.  Journal of Physics A Mathematical General, 42(6):065303, February 2009.
  * [5] R. Alicki, M. Horodecki, P. Horodecki, and R. Horodecki.  On thermal stability of topological qubit in Kitaev’s 4D model.  arXiv e-prints, page arXiv:0811.0033, November 2008.
  * [6] Alexei Kitaev.  Anyons in an exactly solved model and beyond.  Annals of Physics, 321(1):2–111, January 2006.
  * [7] Ahmed Almheiri, Xi Dong, and Daniel Harlow.  Bulk Locality and Quantum Error Correction in AdS/CFT.  JHEP, 04:163, 2015.
  * [8] S. W. Hawking and Don N. Page.  Thermodynamics of black holes in anti-de Sitter space.  Communications in Mathematical Physics, 87(4):577 – 588, 1982.
  * [9] Edward Witten.  Anti-de Sitter space, thermal phase transition, and confinement in gauge theories.  Adv. Theor. Math. Phys., 2:505–532, 1998.
  * [10] Ning Bao, ChunJun Cao, and Guanyu Zhu.  Deconfinement and Error Thresholds in Holography.  arXiv e-prints, page arXiv:2202.04710, February 2022.
  * [11] Fernando Pastawski, Beni Yoshida, Daniel Harlow, and John Preskill.  Holographic quantum error-correcting codes: Toy models for the bulk/boundary correspondence.  JHEP, 06:149, 2015.
  * [12] Gong Cheng, ChunJun Cao, and Brian Swingle.  In preparation.
  * [13] Emanuel Knill, Raymond Laflamme, and Lorenza Viola.  Theory of Quantum Error Correction for General Noise.  Phys. Rev. Lett., 84(11):2525–2528, March 2000.
  * [14] Cédric Bény and Ognyan Oreshkov.  General Conditions for Approximate Quantum Error Correction and Near-Optimal Recovery Channels.  Phys. Rev. Lett., 104(12):120501, March 2010.
  * [15] Eric Mintun, Joseph Polchinski, and Vladimir Rosenhaus.  Bulk-Boundary Duality, Gauge Invariance, and Quantum Error Corrections.  Phys. Rev. Lett., 115(15):151601, 2015.
  * [16] Alexey Milekhin.  Quantum error correction and large N.  SciPost Physics, 11(5):094, November 2021.
  * [17] Eric Dennis, Alexei Kitaev, Andrew Landahl, and John Preskill.  Topological quantum memory.  Journal of Mathematical Physics, 43(9):4452–4505, September 2002\.
  * [18] Dave Bacon.  Operator quantum error-correcting subsystems for self-correcting quantum memories.  Phys. Rev. A, 73:012340, Jan 2006.
  * [19] Benjamin J. Brown, Daniel Loss, Jiannis K. Pachos, Chris N. Self, and James R. Wootton.  Quantum memories at finite temperature.  Reviews of Modern Physics, 88(4):045005, October 2016.
  * [20] Olivier Landon-Cardinal, Beni Yoshida, David Poulin, and John Preskill.  Perturbative instability of quantum memory based on effective long-range interactions.  Phys. Rev. A, 91:032303, Mar 2015.
  * [21] Courtney G. Brell.  A proposal for self-correcting stabilizer quantum memories in 3 dimensions (or slightly less).  New Journal of Physics, 18(1):013050, January 2016.
  * [22] Keiichiro Furuya, Nima Lashkari, and Mudassir Moosa.  Renormalization group and approximate error correction.  12 2021.
  * [23] Robert J. Harris, Nathan A. McMahon, Gavin K. Brennen, and Thomas M. Stace.  Calderbank-Shor-Steane holographic quantum error-correcting codes.  Phys. Rev. A, 98(5):052301, 2018.
  * [24] ChunJun Cao and Brad Lackey.  Approximate Bacon-Shor Code and Holography.  JHEP, 05:127, 2021.
  * [25] Kfir Dolev, Vladimir Calvera, Samuel S. Cree, and Dominic J. Williamson.  Gauging the bulk: generalized gauging maps and holographic codes.  JHEP, 05:158, 2022.
  * [26] Ofer Aharony, Joseph Marsano, Shiraz Minwalla, Kyriakos Papadodimas, and Mark Van Raamsdonk.  The Hagedorn - deconfinement phase transition in weakly coupled large N gauge theories.  Adv. Theor. Math. Phys., 8:603–696, 2004.
  * [27] Nissan Itzhaki and John McGreevy.  The Large N harmonic oscillator as a string theory.  Phys. Rev. D, 71:025003, 2005.
  * [28] David Berenstein.  A toy model for the AdS/CFT correspondence.  Journal of High Energy Physics, 2004(07):018–018, jul 2004.
  * [29] Igor R. Klebanov.  String theory in two-dimensions.  In Spring School on String Theory and Quantum Gravity (to be followed by Workshop), 7 1991.
  * [30] Clifford V. Johnson.  Consistency Conditions for Non-Perturbative Completions of JT Gravity.  12 2021.
  * [31] Stefano Chesi, Beat Röthlisberger, and Daniel Loss.  Self-correcting quantum memory in a thermal environment.  Phys. Rev. A, 82:022305, Aug 2010.
  * [32] Daniel Gottesman, Alexei Kitaev, and John Preskill.  Encoding a qubit in an oscillator.  Physical Review A, 64(1):012310, July 2001.

##  Appendix A Proof of Theorem 1

In this section, we will prove that the mutual information decays polynomially
in $t$ when coupled to a thermal bath of suitably low temperature, as long as
the conditions in Theorem 1 are satisfied. We first set up some useful
notations for later use.

Let $\cal{E}$ be the set of operators that generate the errors. For all
operators $E_{a}\in\mathcal{E}$ satisfying the aKL equation, we can always
diagonalize the matrix $f_{ab}$ by a basis transformation inside
$\mathcal{E}$, such that
$\langle\tilde{i}|E_{b}^{\dagger}E_{a}|\tilde{i}\rangle=f_{a}\delta_{ab}+O(\frac{1}{N^{2}})$.
The matrix elements here are taken with respect to the $d$ dimensional code
words. Then the aKL equation can be written as

$\langle\tilde{j}|E_{b}^{\dagger}E_{a}|\tilde{i}\rangle=f_{a}\delta_{ab}\delta_{ij}+\frac{g_{ab}^{ij}}{N^{2}},$
(92)

where $f_{a}\geq 0$ is some $O(1)$ number that depends on the error indices.

Now we define two classes of errors according to whether $f_{a}$ is 0 or not.
The first class of errors, labeled by $E_{\alpha}$ and $E_{\beta}$, satisfies
$f_{\alpha},f_{\beta}\neq 0$. The second class of errors, labeled by
$E_{\sigma}$ and $E_{\nu}$, satisfies $f_{\sigma}=f_{\nu}=0$. Note that we
will label the errors with Greek indices, $\alpha,\beta$ for the first class,
$\sigma,\nu$ for the second class, when we distinguish their classes. When we
do not need to treat them in different classes, we use Latin indices like
$a,b$ to label these errors. For the first class of errors, we can define an
orthonormal basis $\\{|\mu_{\alpha}^{i}\rangle\\}$ by acting error operators
$E_{\alpha}$ on the logical states $|\tilde{i}\rangle$ such that

$\begin{split}&E_{\alpha}|\tilde{i}\rangle=\sqrt{f_{\alpha}}|\mu_{\alpha}^{i}\rangle+\sum_{\beta}\frac{c_{\alpha\beta}}{2\sqrt{f_{\beta}}N^{2}}|\mu_{\beta}^{i}\rangle+\sum_{\beta}\frac{{\Delta}_{\alpha\beta}^{i}}{2\sqrt{f_{\beta}}N^{2}}|\mu_{\beta}^{i}\rangle+\sum_{\beta,j\neq
i}\frac{h_{\alpha\beta}^{ij}}{2\sqrt{f_{\beta}}N^{2}}|\mu_{\beta}^{j}\rangle.\\\
\end{split}$ (93)

The states $|\mu_{\alpha}^{i}\rangle$’s can be read off from the $O(1)$ term
in above equations. All the other coefficients can be read off from the aKL
equations as,

$\begin{split}&c_{\alpha\beta}=\frac{1}{d}\sum_{i=1}^{d}g_{\alpha\beta}^{ii}\\\
&\Delta^{i}_{\alpha\beta}=g_{\alpha\beta}^{ii}-c_{\alpha\beta}\\\
&h_{\alpha\beta}^{ij}=\begin{cases}g_{\alpha\beta}^{ij},&\text{if $i\neq j$
}\\\ 0,&\text{if $i=j$}\end{cases}.\end{split}$ (94)

One can see that there is a clear separation between the phase errors that are
captured by $\Delta^{i}_{\alpha\beta}$, and the bit flip errors that are
captured by $h_{\alpha\beta}^{ij}$. So far we have discussed the first class
of errors with non-zero $O(1)$ norm $\sqrt{f_{\alpha}}$. However, the second
class of error operators have vanishing $O(1)$ norm and let us label them by
$E_{\sigma}$ and $E_{\nu}$. so $\langle
E_{\sigma}^{\dagger}E_{\sigma}\rangle=O(\frac{1}{N^{2}})$. Generally, the aKL
equations for these operators take the form of

$\langle\tilde{j}|E_{\sigma}^{\dagger}E_{\nu}|\tilde{i}\rangle=\frac{g_{\sigma\nu}^{ij}}{N^{2}}+O(\frac{1}{N^{3}}).$
(95)

These operators in the second class of errors map the logical states to states
with norms proportional to $O(\frac{1}{N})$. More generally, they are
operators that satisfy the following equation:

$E_{\sigma}|\tilde{i}\rangle=\sum_{\nu}\frac{s_{\sigma\nu}^{i}}{N}|\mu_{\nu}^{i}\rangle,$
(96)

for some $O(1)$ coefficients $s$’s. The states $|\mu_{\nu}^{i}\rangle$’s are
another set of normalized basis, orthogonal to the states
$|\mu_{\alpha}^{i}\rangle$’s up to $O(\frac{1}{N})$ overlap. The errors that
satisfy

$\langle\tilde{i}|E_{\sigma}|\tilde{j\rangle}=\frac{e^{ij}_{\sigma}}{N}+O(\frac{1}{N^{3}}),$
(97)

also belong to the second class of errors, since they are consistent with Eq.
(96).

For some calculations, we do not need to differentiate the two types of
errors. In that case, we use Latin indices and simply label all of them by
$E_{a}$. Let us denote the KL matrix element by
$\langle\tilde{j}|E_{b}^{\dagger}E_{a}|\tilde{i}\rangle=G_{ab}^{ij}$; one can
write down a general relation

$\begin{split}&E_{a}|\tilde{i}\rangle=\sum_{A}u_{aA}^{i}|\mu_{A}\rangle,\end{split}$
(98)

where $\\{|\mu_{A}\rangle\\}$ is a set of orthonormal states, and $u^{i}_{aA}$
are coefficients that satisfy

$\sum_{A}u_{aA}^{i}u_{bA}^{j}=G_{ab}^{ij}.$ (99)

###  A.1 Early time approximation

Recall that we model our open system dynamics as a Markovian process. At early
time, we approximate the evolution by only including the first order
perturbation of errors.141414The coupling $\lambda_{a}(\epsilon_{a})$ for each
error operator $E_{a}$ is defined as $\lambda_{a}\gamma(\epsilon_{a})$,, with
$\gamma(\epsilon_{a})$ defined in Eq.22.

$\rho(t)\approx\rho_{0}-\sum_{a}\frac{\lambda_{a}(\epsilon_{a})t}{2}\\{E_{a}^{\dagger}E_{a},\rho_{0}\\}+\sum_{a}\lambda_{a}(\epsilon_{a})tE_{a}\rho_{0}E_{a}^{\dagger}+\cdots,$
(100)

where the sum over $E_{a}$ includes both the first and second class of errors,
and $\rho_{0}$ is the maximally entangled state between reference $R$ and
system $S$ at $t=0$,

$\rho_{0}=\frac{1}{d}\sum_{ij}|\tilde{i}i\rangle_{RS}\langle\tilde{j}j|.$
(101)

We proceed to compute the mutual information
$I(R:S)=S(\rho^{R})+S(\rho^{S})-S(\rho)$ at early times. Since the non-trivial
unitary dynamics that induces the above Markovian evolution only occurs
between the system and the bath, the entropy of $R$ should remain invariant
over time, i.e. $S(\rho^{R})=\ln d$. Therefore, it suffices to calculate the
entropies of $S$ and the $RS$ joint system.

Generally, the density matrix is block-diagonal. We can generate each block
matrix using the operators $E_{a}$, with $a\in\mathscr{B}_{n}$, such that for
any $a,b\in\mathscr{B}_{n}$, $G_{ab}\neq 0$. Here $\mathscr{B}_{n}$ is an
index set that captures the error operators used to generate a particular
block matrix. Sometimes we will also refer to this set as a block, since it
defines the corresponding block matrix. We write the density matrix associated
to the block $\mathscr{B}_{n}$ as,

$\begin{split}\rho_{n}=&\sum_{a\in\mathscr{B}_{n}}E_{a}\rho_{0}E^{\dagger}_{a}\\\
=&\frac{1}{d}\sum_{aABij}u_{aA}^{i}u_{aB}^{j*}|\mu_{A}\rangle_{S}\langle\mu_{B}|\otimes|i\rangle_{R}\langle
j|.\end{split}$ (102)

The eigenvalues of matrix $\frac{1}{d}\sum_{a}u_{aA}^{i}u_{aB}^{j*}$ is the
same as that of the matrix
$\frac{1}{d}\sum_{iA}u_{aA}^{i}u_{bA}^{i*}=\frac{1}{d}\sum_{i}G_{ab}^{ii}$.
After tracing out the reference $R$, the density matrix of system becomes

$\rho_{n}^{S}=\frac{1}{d}\sum_{aiAB}u_{aA}^{i}u_{aB}^{i*}|\mu_{A}\rangle_{S}\langle\mu_{B}|.$
(103)

The eigenvalues of matrix $\frac{1}{d}\sum_{ai}u_{aA}^{i}u_{aB}^{i*}$ is the
same as that of the matrix
$\frac{1}{d}\sum_{A}u_{aA}^{i}u_{bA}^{j*}=\frac{1}{d}G_{ab}^{ij}$.

Therefore, to evaluate the entropy, we need to calculate the eigenvalues of
the matrix $\frac{1}{d}G_{ab}^{ij}$. Generally, inside a block, we have both
the first and the second class of errors. Their overlap $G_{\alpha\sigma}$
scales as $\frac{1}{N^{2}}$ by the aKL condition. However, we can always
diagonalize the matrix $G$, where the diagonal form $\tilde{G}$ has vanishing
matrix elements for the overlaps $\tilde{G}_{\alpha\sigma}=0$ and
$\tilde{G}_{\alpha\beta}-G_{\alpha\beta}=O(\frac{1}{N^{4}}),\tilde{G}_{\sigma\nu}-G_{\sigma\nu}=O(\frac{1}{N^{4}})$
such that the correction is at subleading order. Therefore, in the following,
we can neglect the overlap between the first and the second class of errors,
and always treat them as different blocks.

For our following discussion, we use $\rho_{n}^{(1)},\rho_{n}^{(2)}$ to denote
the state created purely by the first and the second class of errors,
respectively. The entropy for $\rho_{n}^{(1)}$ is,

$\begin{split}S(\rho_{n}^{(1)})&=-\operatorname{Tr}(\rho_{n}^{(1)}\ln\rho_{n}^{(1)})\\\
&=-\operatorname{Tr}[(f+\frac{c}{N^{2}})\ln(f+\frac{c}{N^{2}})]\\\
&=-\operatorname{Tr}(f\ln f)-\frac{1}{N^{2}}\operatorname{Tr}(c\ln
f)-\frac{1}{N^{2}}\operatorname{Tr}(c)-\frac{1}{2N^{4}}\operatorname{Tr}(cf^{-1}c)+O(\frac{1}{N^{6}}).\end{split}$
(104)

In the second step we used that the eigenvalues of the density matrix are the
same as those of the KL overlap matrix $\frac{1}{d}\sum_{i}G^{ii}$. $f$ is the
diagonal matrix $f_{\alpha}\delta_{\alpha\beta}$ and the matrix elements of
$c$ are given by Eq. (94). For the reduced density matrix of system $S$,
$\rho^{S}_{n}=\operatorname{Tr}_{R}(\rho_{n})$, its entropy is

$\begin{split}S(\rho_{n}^{(1)S})=&-\operatorname{Tr}[\frac{1}{d}(f\otimes
I+\frac{g}{N^{2}})\ln(f\otimes
I+\frac{g}{N^{2}})]+\frac{1}{d}\operatorname{Tr}(f\otimes
I+\frac{g}{N^{2}})\ln{d}\\\
=&\operatorname{Tr}(\rho_{n})\ln{d}-\frac{1}{d}\operatorname{Tr}[(f\otimes
I)\ln(f\otimes I)]-\frac{1}{dN^{2}}\operatorname{Tr}[g\ln(f\otimes
I)]-\frac{1}{dN^{2}}\operatorname{Tr}(g)\\\
&-\frac{1}{2dN^{4}}\operatorname{Tr}[g(f^{-1}\otimes
I)g]+O(\frac{1}{N^{6}}).\end{split}$ (105)

where $f$ is the matrix with element being $f_{\alpha\beta}$, and $I$ is the
$d\times d$ identity matrix. After subtracting $S(\rho_{n}^{(1)})$ from
$S(\rho_{n}^{(1)S})$, we finally obtain

$\begin{split}S(\rho_{n}^{(1)S})-S(\rho_{n}^{(1)})&=\operatorname{Tr}(\rho_{n}^{(1)})\ln{d}-\frac{1}{2dN^{4}}\operatorname{Tr}[g(f^{-1}\otimes
I)g]+\frac{1}{2N^{4}}\operatorname{Tr}(cf^{-1}c)+O(\frac{1}{N^{6}})\\\
&=\operatorname{Tr}(\rho_{n}^{(1)})\ln
d-\frac{1}{2dN^{4}}(\sum_{i}\operatorname{Tr}[\Delta^{i}f^{-1}\Delta^{i}]+\sum_{ij}\operatorname{Tr}[h_{ij}f^{-1}h_{ij}^{T}])+O(\frac{1}{N^{6}}),\end{split}$
(106)

where matrices $\Delta^{i}$ and $h_{ij}:=h^{ij}$ are defined in Eq. (94). We
see that the decay rate of mutual information is suppressed by
$\frac{1}{N^{4}}$. Repeating the same subtractions but for the second class of
errors, the difference in entropy is

$\begin{split}S(\rho_{n}^{(2)S})-S(\rho_{n}^{(2)})=&\frac{1}{N^{2}}[S(\frac{1}{d}g^{ij}_{\sigma\nu})-S({\frac{1}{d}\sum_{i}g^{ii}_{\sigma\nu}})]|_{\sigma,\nu\in\mathscr{B}_{n}}+O(\frac{1}{N^{4}})\\\
=&\operatorname{Tr}(\rho_{n}^{(2)})\ln{d}-\frac{1}{N^{2}}[S(\frac{1}{d}{\sum_{i}g^{ii}_{\sigma\nu}})-\frac{1}{d}S(g^{ij}_{\sigma\nu})]|_{\sigma,\nu\in\mathscr{B}_{n}^{\prime}}+O(\frac{1}{N^{4}}),\end{split}$
(107)

where $S(g_{\sigma\nu}^{ij})$ means that we are calculating entropy of the
matrix whose elements are $g_{\sigma i,\nu j}$. On the other hand, in
$S(\frac{1}{d}\sum_{i}g_{\sigma\nu}^{ii})$, we compute the entropy of the
matrix with elements proportional to $\frac{1}{d}\sum_{i}g^{ii}_{\sigma\nu}$.
The decay rate is only suppressed in $\frac{1}{N^{2}}$. The entropy difference
in the second term depends on the matrix elements $g_{\sigma\nu}^{ij}$ in a
complicated way. Instead of calculating it explicitly, we just provide an
upper bound for this expression. Suppose that in each block $\mathscr{B}_{n}$,
the matrix $g_{\sigma\nu}^{ij}$ has a rank of $r(n)$ with its elements bounded
by a function $f(n)$, then the second term is bounded by
$\frac{r(n)f(n)\log{d}}{N^{2}}$.

For $n\geq 1$, define $\alpha_{n}=\lambda_{n}t$. The early time density matrix
can be written as151515To keep the notation tractable, we absorbed any
$a$-dependence of the coupling $\lambda_{a}$ into the operators $E_{a}$ such
that the couplings $\lambda$ only retain the $n$ dependence. To be more
precise, one can redefine the couplings and errors such that
$\lambda(\epsilon_{a})$ only depends on the energy of each error operator
$\epsilon_{a}$. Because $\epsilon_{a}=\epsilon_{n}$ for all
$a\in\mathscr{B}_{n}$, we can rewrite $\lambda(\epsilon_{a})$ as
$\lambda_{n}$.

$\rho(t)=\bigoplus_{n=0}^{\infty}(\alpha_{n}\rho_{n}).$ (108)

Applying the relation $\sum_{n}\alpha_{n}\operatorname{Tr}(\rho_{n})=1$, one
find that the difference of total entropy is

$\begin{split}S(\rho^{S}(t))-S(\rho(t))=&\ln{d}-\sum_{n=1}^{\infty}\frac{\lambda_{n}t}{N^{2}}[S({\frac{1}{d}\sum_{i}g^{ii}_{\sigma\nu}})-\frac{1}{d}S(g^{ij}_{\sigma\nu})]|_{\sigma,\nu\in\mathscr{B}_{n}^{\prime}}\\\
&-\sum_{n=1}^{\infty}\frac{\lambda_{n}t}{2dN^{4}}\sum_{\alpha\beta\in\mathscr{B}_{n}}[\sum_{i}\Delta^{i2}_{\alpha\beta}+\sum_{i\neq
j}(h_{\alpha\beta}^{ij})^{2}]+O(\frac{1}{N^{6}}).\end{split}$ (109)

For early time, i.e., $\lambda_{n}t\ll 1$, the dominant contribution to
information decay is from the second class of errors.

###  A.2 Finite time

In this section, we look at the finite time behavior of density matrix up to
$O(\frac{1}{N^{2}})$ contributions. We will see that, although the early-time
effect of the first class of errors is sub-dominant, it becomes more relevant
at later times since it grows quadratically in time. Recall that the overlap
between the different classes of errors are suppressed, we can split the
Lindbladian operator into two parts associated with the first class and the
second class of errors,

$\rho(t)=e^{\mathcal{L}t}\rho(0)=e^{\mathcal{L}^{(1)}t+\mathcal{L}^{(2)}t}\rho(0)+O(\frac{1}{N^{4}}).$
(110)

From Eq. (96), the second class of Lindblad operator creates a term of
$O(\frac{1}{N^{2}})$ each time it acts on the density matrix. Therefore, the
state can be approximated by

$\rho(t)=\rho^{(1)}(t)+t\mathcal{L}^{(2)}[\rho^{(1)}(t)]+O(\frac{1}{N^{4}}),$
(111)

with $\rho^{(1)}(t)=e^{\mathcal{L}^{(1)}t}\rho(0)$ being the state evolved
under the first class of errors. Here we only keep the commuting contributions
from $\mathcal{L}^{(1)}$ and $\mathcal{L}^{(2)}$, since their commutator only
contributes at higher orders of $\frac{1}{N}$.

To simplify the discussion, we will look at an approximation, which we claim
is sufficient to capture the behavior of information decay in our discussion.
Let us first look at the case with only one type of error $E$ in the first
class,

$\dot{\rho}(t)=\lambda E\rho E^{\dagger}+\lambda E^{\dagger}\rho
E-\frac{\lambda}{2}\\{E^{\dagger}E,\rho\\}-\frac{\lambda}{2}\\{EE^{\dagger},\rho\\}.$
(112)

For low energy operators $E$, we can make the approximation where the coupling
of $E$ and its conjugate has the same strength $\lambda$, which in general
should be different due to the thermal factor. In this large $N$ model, we
observe that $E^{\dagger}E$ acts trivially in the leading order,

$\begin{split}E^{\dagger}EE^{p}\rho&=(p+1)E^{p}\rho+O(\frac{1}{N^{2}})\\\
EE^{\dagger}E^{p}\rho&=pE^{p}\rho+O(\frac{1}{N^{2}}).\end{split}$ (113)

where $p$ is some integer power.

The order $\frac{1}{N^{2}}$ correctons can be captured by an operator $Z$,
which acts on the code state in the following way:

$Z|\mu^{i}_{\alpha}\rangle=(1+\frac{g^{ii}}{2N^{2}})|\mu^{i}_{\alpha}\rangle.$
(114)

It is not hard to find the solution up to order $\frac{1}{N^{2}}$,

$\rho(t)=\sum_{k=0}^{\infty}\frac{(\lambda t)^{k}}{(1+\lambda
Zt)^{k+1}}\frac{E^{k}\rho_{0}E^{\dagger k}}{k!}+O(\frac{1}{N^{4}}).$ (115)

One might generalize this equation to include all the first class error
operators $E_{\alpha}$. We denote the density matrix dressed by the first
class of errors as $\rho^{(1)}(t)$, which is

$\rho^{(1)}(t)=\sum_{\\{k_{\alpha}\\}}w(\\{k_{\alpha}\\})\prod_{\alpha}\frac{E_{\alpha}^{k_{\alpha}}\rho(0)E_{\alpha}^{\dagger
k_{\alpha}}}{k_{\alpha}!}+O(\frac{1}{N^{4}}),$ (116)

where the weight function can be upper bounded by,

$\begin{split}&w(\\{k_{\alpha}\\})=\prod_{\alpha}w(k_{\alpha}),\\\
&\mathrm{where}~{}w(k_{\alpha})\leq\frac{(\lambda_{\alpha}t)^{k_{\alpha}}}{(1+\lambda_{\alpha}Z_{\alpha}t)^{k_{\alpha}+1}},\
\ \text{for $k_{\alpha}\geq 1$}.\end{split}$ (117)

The bound is only saturated if we could neglect the difference between the
coupling constants of $E_{\alpha}$ and $E_{\alpha}^{\dagger}$. While said
difference is negligible for low energy errors, it is not for high energy
ones. In fact, the coupling constant goes down exponentially with the energy
$\epsilon_{\alpha}$ of the operator $E_{\alpha}$,

$\lambda_{\alpha}\propto\gamma(\epsilon_{\alpha})\sim
e^{-\beta\epsilon_{\alpha}},$ (118)

for $\epsilon_{\alpha}>0$. However, the coupling constant
$\lambda_{\alpha}^{\prime}$ for $E_{\alpha}^{\dagger}$ barely depends on
energy, since $\gamma(\epsilon_{\alpha})\sim 1$, for $\epsilon_{\alpha}<0$. We
only bound the weight function for $k_{\alpha}\geq 1$, since the
$k_{\alpha}=0$ component does not contribute to the mutual information decay,
and it is determined by the constraint
$\sum_{k_{\alpha}=0}^{\infty}w(k_{\alpha})=1$ of the weight function.

Given the density matrix, we calculate the entropy difference by the equation,

$S(\rho^{(1)S})-S(\rho^{(1)})=\operatorname{Tr}(\rho^{(1)})\ln{2}-\frac{1}{2dN^{4}}\sum_{\mathcal{A}\mathcal{B}}w_{\mathcal{A}}(\sum_{i}(\Delta_{\mathcal{A}\mathcal{B}}^{i})^{2}+\sum_{i\neq
j}(h_{\mathcal{A}\mathcal{B}}^{ij})^{2})+O(\frac{1}{N^{6}}),$ (119)

where we use $\mathcal{A}$ and $\mathcal{B}$ as a shorthand notation for the
product of first class operators, so that
$E_{\mathcal{A}}:=\prod_{\alpha_{i}\in\mathcal{A}}E_{\alpha_{i}}$.
$w_{\mathcal{A}}$ is a number obtained by setting all the $Z_{\alpha}$’s to
identity in Eq. (117). The $d$-function and $h$-function scale linearly with
$k_{n}$, and are proportional to some power $q$ of $n_{\alpha}$, with
$n_{\alpha}=\frac{\epsilon_{\alpha}}{\Omega}$ being the energy level of
$E_{\alpha}$, where $\Omega$ has dimension of energy. For example,
$\Omega=\omega$ for our oscillator model. Thus we can upper bound the
information decay rate by (up to $\frac{1}{N^{4}}$),

$\begin{split}K^{(1)}(t)&=\frac{1}{N^{4}}\sum_{k_{\alpha}=1}^{\infty}w(\\{k_{\alpha}\\})(\sum_{\alpha}n_{\alpha}^{q}k_{\alpha})^{2}\\\
&\leq\frac{1}{N^{4}}[2(\sum_{\alpha}\lambda_{\alpha}n_{\alpha}^{q}t)^{2}+\sum_{\alpha}\lambda_{\alpha}n_{\alpha}^{2q}t].\end{split}$
(120)

Note that the coefficient $\lambda_{\alpha}$ scales as
$e^{-\beta\epsilon_{\alpha}}$, and the density of errors with energy
$\epsilon$ is bounded by $e^{\mu\epsilon}$ by our sparsity condition. By the
same argument we have for Eq.(58) in Section 3.2, one can show that the
quantity inside the square bracket does not scale as $N$ as long as
$T<\frac{1}{\mu}$.

The decay rate due to the first class of errors is suppressed by $N^{4}$ and
grows quadratically with time for $\lambda t\gg 1$. Similar arguments also
apply for higher orders of $\frac{1}{N^{2}}$. In general, the decay rate at
order $\frac{1}{N^{2m}}$ is proportional to $\Delta^{m}$ and $h^{m}$ (defined
in Eq. (94)) of the corresponding error operators, which scales as $k^{m}$ for
the error operator $E_{\alpha}^{k}$. Therefore, the general decay rate takes
the form of $\frac{t^{m}}{N^{2m}}$ at each order with $m\geq 2$.

The second class of errors of the form $E_{\nu}^{k}$ is suppressed by
$\frac{t^{k}}{N^{2k}}$ in the Lindblad equation, where the leading order term
is linear in $t$,

$\begin{split}\rho^{(2)}(t)&:=t\mathcal{L}^{(2)}[\rho^{(1)}(t)]\\\
&=\sum_{\nu}\lambda_{\nu}t[E_{\nu}\rho^{(1)}(t)E_{\nu}^{\dagger}-\frac{1}{2}\\{\rho^{(1)}(t),E_{\nu}^{\dagger}E_{\nu}\\}].\end{split}$
(121)

Recall that $\rho^{(1)}(t)$ is the density matrix induced by the first class
of errors, explicitly given by Eq. (116). For the sake of clarity, let us
rewrite $\rho^{(1)}(t)$ as

$\rho^{(1)}(t)=\sum_{\mathcal{A}}w_{\mathcal{A}}(t)E_{\mathcal{A}}\rho_{0}E_{\mathcal{A}}^{\dagger}.$
(122)

Like our discussion in Appendix A.1, we convert the entropy of the state to
the entropy of the matrix $g_{\sigma\mathcal{A},\nu\mathcal{B}}^{ij}$, which
can be read off from

$\langle\tilde{j}|E^{\dagger}_{\mathcal{B}}E^{\dagger}_{\nu}E_{\sigma}E_{\mathcal{A}}|\tilde{i}\rangle=\frac{g_{\mathcal{A}\sigma,\mathcal{B}\nu}^{ij}}{N^{2}}+O(\frac{1}{N^{4}}).$
(123)

Being the product of a second class error and a first class error, this
composite operator belongs to the second class of errors, but it has
$O(\frac{1}{N})$ overlap with some operators in the first class,

$\langle\tilde{j}|E_{\mathcal{B}}^{\dagger}E_{\sigma}E_{\mathcal{A}}|\tilde{i}\rangle=\frac{\delta_{ij}g_{\sigma\mathcal{A},\mathcal{B}}}{N}+O(\frac{1}{N^{3}}).$
(124)

Due to the factorization property,
$\langle\tilde{j}|E^{\dagger}_{\mathcal{B}}E^{\dagger}_{\nu}E_{\sigma}E_{\mathcal{A}}|\tilde{i}\rangle=\langle\tilde{j}|E_{\nu}^{\dagger}E_{\sigma}|\tilde{i}\rangle\langle
0|E_{\mathcal{B}}^{\dagger}E_{\mathcal{A}}|0\rangle+\langle\tilde{j}|\tilde{i}\rangle\langle
0|E_{\mathcal{B}}^{\dagger}E_{\nu}^{\dagger}E_{\sigma}E_{\mathcal{A}}|0\rangle$,
the matrix element can be written as sum of two terms,

$g^{ij}_{\mathcal{A}\sigma,\mathcal{B}\nu}=g^{ij}_{\sigma\nu}f_{\mathcal{A}\mathcal{B}}+\delta^{ij}g_{\sigma\mathcal{A},\nu\mathcal{B}}.$
(125)

where $f_{\mathcal{A}\mathcal{B}}:=\langle
0|E_{\mathcal{B}}^{\dagger}E_{\mathcal{A}}|0\rangle$, and $|0\rangle$ is the
vacuum state of our system. The matrix $g^{ij}_{\sigma\nu}$ is identical to
the one in Eq. (95). $g_{\sigma\mathcal{A},\nu\mathcal{B}}$ is positive
definite where each element is equal to the four point function $\langle
0|E_{\mathcal{B}}^{\dagger}E_{\nu}^{\dagger}E_{\sigma}E_{\mathcal{A}}|0\rangle$,
modulo the $\frac{1}{N^{2}}$ factor. The overlap between the composite
operators and some first class operators (Eq. (124)) can be removed by a basis
transformation, which also cancels the second term in Eq. (125). Therefore,
similar to Eq. (107), we conclude that the entropy difference is

$\begin{split}S(\rho^{(2)S}(t))-S(\rho^{(2)}(t))=&\frac{1}{N^{2}}[S(\frac{1}{d}g^{ij}_{\sigma\nu}(t)f_{\mathcal{A}\mathcal{B}}(t))-S(\frac{1}{d}\sum_{i}g^{ii}_{\sigma\nu}(t)f_{\mathcal{A}\mathcal{B}}(t))]\\\
=&\frac{\operatorname{Tr}(f_{\mathcal{A}\mathcal{B}}(t))}{N^{2}}[S(\frac{1}{d}g_{\sigma\nu}^{ij}(t))-S(\frac{1}{d}\sum_{i}g_{\sigma\nu}^{ii}(t))]\\\
=&\sum_{n}\frac{\lambda_{n}t}{N^{2}}[S(\frac{1}{d}g_{\sigma\nu}^{ij})-S(\frac{1}{d}\sum_{i}g_{\sigma\nu}^{ii})]|_{\sigma,\nu\in\mathscr{B}_{n}^{\prime}},\end{split}$
(126)

where $S(\frac{1}{d}g_{\sigma\nu}^{ij}(t)f_{\mathcal{A}\mathcal{B}}(t))$
denotes the entropy of the matrix whose element is
$\frac{1}{d}g_{\sigma\mathcal{A}i,\nu\mathcal{B}j}(t)=\frac{1}{d}g_{\sigma\nu}^{ij}(t)f_{\mathcal{A}\mathcal{B}}(t)$.
We have absorbed the factor $\lambda_{\nu}t$ into the matrix
$g_{\sigma\nu}^{ij}(t):=\sqrt{\lambda_{\nu}\lambda_{\sigma}}tg_{\sigma\nu}^{ij}$,
and the weight function $w_{\mathcal{A}}(t)$ into the function
$f_{\mathcal{A}\mathcal{B}}(t):=\sqrt{w_{\mathcal{A}}(t)w_{\mathcal{B}}(t)}f_{\mathcal{A}\mathcal{B}}$.
By definition,
$\operatorname{Tr}(f_{\mathcal{A}\mathcal{B}})=\sum_{\mathcal{A}}w_{\mathcal{A}}(t)\langle
0|E_{\mathcal{A}}^{\dagger}E_{\mathcal{A}}|0\rangle=\operatorname{Tr}(\rho^{(1)}(t))+O(\frac{1}{N^{2}})=1+O(\frac{1}{N^{2}})$.
In the last step, we again utilized the block diagonal form of the matrix
$g_{\sigma\nu}^{ij}$ to write the result as a sum of each block’s
contribution, with blocks being labeled by $\mathscr{B}_{n}^{\prime}$, such
that the matrix elements are non-zero only for $\sigma$ and $\nu$ in the same
block. Note that the entropy difference is the same as in Eq. (107). We can
bound the decay rate induced by the second class of errors with the following
function at leading order of $\frac{\lambda t}{N^{2}}$,

$K^{(2)}(t)\leq\sum_{n}\frac{\lambda_{n}t}{N^{2}}f(n)r(n)\log{d},$ (127)