[1]School of Physics, Georgia Institute of Technology, Atlanta, GA, USA.
[2]Department of Physics, The University of Texas at Austin, Austin, TX, USA.

# Matrix product states and the decay of quantum conditional mutual
information

Pavel Svetlichnyy  Shivan Mittal  T.A.B. Kennedy

###### Abstract

A uniform matrix product state defined on a tripartite system of spins,
denoted by $ABC,$ is shown to be an approximate quantum Markov chain when the
size of subsystem $B,$ denoted $|B|,$ is large enough. The quantum conditional
mutual information (QCMI) is investigated and proved to be bounded by an
exponentially decaying function of $|B|$ in the same limit. The properties of
the bounding function are derived by a new approach, with a corresponding
improved value given for its asymptotic decay rate. Numerical investigations
of the decay of QCMI are reported for a collection of matrix product states
generated by selecting the defining isometry with respect to Haar measure.

##  1 Introduction

Consider the collection of subsets of a set $X,$ say the configuration space
of a physical system, and a positive-valued set function $f$ on $X$. Then $f$
is said to be strongly subadditive if $f(U)+f(V)\geq f(U\cap V)+f(U\cup V),$
for any subsets $U$ and $V$ of $X$. Suppose $A$ , $B$ and $C$ are mutually
disjoint subsets and let $U=A\cup B=:AB$ and $V=B\cup C=:BC,$ so that $U\cup
V=ABC$ and $U\cap V=B$ [Wil17]. Then for $f$ strongly subadditive

$f(AB)+f(BC)\geq f(B)+f(ABC).$

We focus attention on a one dimensional, homogeneous spin chain that is
subdivided into a tripartite system with regions $A$ and $C$ separated by
“buffer region” $B.$ For any such configuration space region, denoted by $R$
with Hilbert space $\mathcal{H}_{R},$ the density operator is denoted by
$\rho_{R}\in\mathcal{B}(\mathcal{H}_{R}).$ Then the quantum entropy
$S(\rho_{R})=-\mbox{tr}\left(\rho_{R}\ln\rho_{R}\right),$ is such a strongly
subadditive function, as conjectured by Lanford and Robinson [LR68] and proved
by Lieb and Ruskai [LR73, LR73a].

The strong subadditivity of quantum entropy theorem states that for any
density operator $\rho_{ABC}$ on Hilbert space
$\mathcal{H}_{A}\otimes\mathcal{H}_{B}\otimes\mathcal{H}_{C}$, the quantum
conditional mutual information (QCMI) of subsystems $A$ and $C$ given $B$,
denoted $I(A:C|B),$111The QCMI may be contrasted with the information mutually
shared by regions $A$ and $C$ and described by the quantum mutual information
(QMI) $I(A:C):=S(\rho_{A})+S(\rho_{C})-S(\rho_{AC}).$ is non-negative

$\displaystyle I(A:C|B)$
$\displaystyle\;\;:=S(\rho_{AB})+S(\rho_{CB})-S(\rho_{ABC})-S(\rho_{B})\geq
0.$

The states $\rho_{ABC}$ for which equality holds are the quantum Markov chains
[Hay+04], and these have the property that they may be recovered exactly from
a reduced density operator; that is, there exists a completely positive map
which induces the transformation $\rho_{AB}\mapsto\rho_{ABC}$ [Pet86, Pet88,
Hay+04], by acting on the spins in the separating region $B$ only. While exact
recovery is associated with the quantum states for which $I(A:C|B)=0$
[Hay+04], approximate recovery requires investigation of states for which QCMI
may be small $I(A:C|B)<g(\varepsilon),$ where $g$ is a known positive function
and $\varepsilon:=\|\rho_{ABC}-\mathcal{R}_{B\rightarrow BC}(\rho_{AB})\|_{1}$
determines the recovery error for the recovery map $\mathcal{R}_{B\rightarrow
BC}$ with domain the support of the reduced density operator $\rho_{B}$ [FR15,
SFR16, Jun+18, Sut18]. Conversely, there exists a so-called universal recovery
map $\mathcal{R}^{u}_{B\rightarrow BC}$ which guarantees that when $I(A:C|B)$
is small, the recovery error is too [FR15, SFR16, Jun+18, Sut18]

$\|\rho_{ABC}-\mathcal{R}^{u}_{B\rightarrow
BC}(\rho_{AB})\|_{1}\leq\sqrt{(4\ln 2)I(A:C|B)}.$

The recovery of quantum states from a quantum memory that has been locally
corrupted in region $C$, say, is thus bounded by the QCMI $I(A:C|B)$ when the
recovery process implements the universal recovery map
$\mathcal{R}^{u}_{B\rightarrow BC}(\rho_{AB})$. It follows, in the case that
$I(A:C|B)$ is a rapidly decaying function of the size of the buffer region,
the number of spins $|B|,$ that the quantum state may be globally recovered by
experiments on the well-localized interface region $B$ and $C$ [Pet86, Pet88,
Hay+04, FR15, SFR16, Jun+18, Fla+17, KB19, BK19, SM16, HP20].

The many-body ground states of spin lattices described by local gapped
Hamiltonians are well approximated by matrix product states (MPS) in one
dimension [KSZ93, Vid03, Pér+07], and their generalization, the projected
entangled pair states (PEPS) in higher dimensions [VMC08]. The MPS first
effectively appeared as finitely correlated states (FCS) [FNW92, FNW94] and
were defined on infinite one-dimensional spin chains. (More precisely, MPS
correspond to the purely generated FCS (pgFCS), which comprise a subclass of
FCS.) The MPS provide an alternative description of the previously developed
density matrix renormalization group algorithm [Whi92], and are also useful in
time dependent studies of weakly entangled systems. In one dimension, MPS may
be generalized to matrix product density operators (MPDO) [VGC04], which, in
the translationally invariant infinite case correspond to the FCS. The MPDO
may be used to approximate thermal states of local gapped Hamiltonians [VGC04,
KB19, KKB20, KAA21]. Controlled experimental preparation techniques for MPS
have been proposed in the literature [Sch+07, ZCL21].

Given the significance of ground and thermal states of local gapped
Hamiltonians, their quantum-information theoretic properties are of great
interest [SM16, Wol+08, BH15, Gut+20]. As are the properties of their
approximants – the MPS and MPDO, which are often easier to study, given their
more explicit structure [Kim17, CKB21, SK22].

In a recent paper [SK22] we proved that the pgFCS (infinite uniform MPS
(iuMPS)), with reduced density operators $\rho_{ABC}$, are approximate quantum
Markov chains for the one-dimensional configuration space $ABC$ discussed
above. Related to this result, the QCMI $I(A:C|B)$ decays and is bounded by a
(sub)exponentially decaying function in the size of the separating region $B$
[SK22]. The goal of the current paper is to reexamine this problem and report
an improved prediction for the exponential decay rate, differing by a factor
of four from our earlier prediction. We also present numerical experiments in
support of this result.

Our method of estimating the QCMI builds on earlier works [FNW92, FNW94,
Pér+07, BH15] and [SK22], where we construct a quantum Markov chain
$\tilde{\rho}_{ABC}$ approximating $\rho_{ABC}$. Specifically, the method is
based on (i) bounding the QCMI by the trace distance between the exact density
operator $\rho_{ABC}$ and its approximant $\tilde{\rho}_{ABC}$, and (ii)
bounding this trace distance by the distance in diamond norm between the
induced quantum channels $\mathcal{E}^{|B|}$ and $\tilde{\mathcal{E}}^{|B|}$
(our notation is explained below). In our earlier work [SK22] the combination
of Alicki-Fannes inequality [AF04, Win16] and Strinespring Continuity theorem
[Sti55, KSW08] provided (i) a QCMI bound linear in trace distance, and (ii)
the trace distance bounded by the square root of the channel distance. In the
present work the QCMI bound is improved to a quadratic dependence on trace
distance, while the latter is bounded by a linear function of the channel
distance. To achieve these improvements, use of the Alicki-Fannes Inequality
has been replaced by a bound based on the Taylor expansion of a logarithmic
function representation of the density operator, with estimate of the
remainder. The continuity of Stinespring dilation [KSW08] is no longer
employed and instead a construction used in the proof of Choi’s theorem
[Cho75] takes advantage of the nearness of the quantum channels, and provides
a sharper tool in this context. Overall we go from a bound on QCMI which had
square root dependence in the quantum channel distance, to one with quadratic
dependence. This leads to a factor of four change in our predicted decay rate,
a result which appears to be corroborated by our numerical experiments using
tensor network methods.

The rest of the paper is organized as follows. In Section 2 we introduce the
framework for MPS and iuMPS. In Section 3 we present and prove our main
result, Theorem 1. In Section 4 we describe and discuss the methods used in
our numerical experiments, and present numerical results for comparison with
Theorem 1 in Section 5. Section 6 contains concluding remarks. A series of
appendices contains supporting technical information.

$M$$M^{\dagger}$$\sigma$$=$$\sigma$(b)$M$$M^{\dagger}$$=$(a)$M$$M^{\dagger}$$E_{jj^{\prime};ii^{\prime}}$
$=$$j$$i$$j^{\prime}$$i^{\prime}$(c)$M$$M^{\dagger}$$X_{k}$$=$
$\nu_{k}$$X_{k}$(d) Figure 1: Diagrammatic depictions of: (a) the iuMPS in
canonical form, as expressed in (2); (b) the density operator $\sigma$, the
fixed point of the quantum channel $\mathcal{E}$, defined in (4); (c) the
transfer matrix $E$ defined in (5); (d) the (right) eigenvalue equation (6)
for $E$.

##  2 Theoretical formulation

A pure state vector for a system of $n$ spins, each with associated Hilbert
space $\mathcal{H}_{s}$ and orthonormal basis
$\\{|s\rangle:s=1,2,..,d_{s}\\}$, may be expressed in the form

$|\Psi\rangle=\sum_{s_{1},\dots,s_{n}=1}^{d_{s}}a_{s_{1}s_{2}\dots
s_{n}}|s_{n}\rangle\otimes\dots|s_{2}\rangle\otimes|s_{1}\rangle,$

where $d_{s}$ is the dimension of $\mathcal{H}_{s}$. For the class of states
called matrix product states (MPS) each coefficient $a_{s_{1}s_{2}\dots
s_{n}}$ has the form [Vid03, Pér+07]

$a_{s_{1}s_{2}\dots s_{n}}=\langle L|M^{s_{n}}[n]\dots
M^{s_{2}}[2]M^{s_{1}}[1]|R\rangle,$

where each $M^{s_{k}}[k]$ belongs to the set
$\mathcal{S}_{k}:=\\{M^{s}[k]\in\mathbb{C}^{d_{M}\times
d_{M}}\;|\;s=1,2\dots,d_{s}\\}$ of complex matrices, indexed by
$k=1,2,\dots,n$, and $|L\rangle,|R\rangle\in\mathbb{C}^{d_{M}}$. We will refer
to the space $\mathcal{H}_{M}:=\mathbb{C}^{d_{M}}$ as virtual space, and to
the space $\mathcal{H}_{s}$ as physical space. For uniform MPS (uMPS) the set
$\mathcal{S}_{k}$ is the same for each site in the chain,
$\mathcal{S}_{k}=\mathcal{S}:=\\{M^{s}\in\mathbb{C}^{d_{M}\times
d_{M}}\;|\;s=1,2\dots,d_{s}\\}$ for $k=1,2,\dots,n$, and uMPS $|\Psi\rangle$
may be expressed as

$|\Psi\rangle=\sum_{s_{1},s_{2},\dots s_{n}=1}^{d_{s}}\langle L|M^{s_{n}}\dots
M^{s_{2}}M^{s_{1}}|R\rangle\;|s_{n}\rangle\otimes\dots\otimes|s_{2}\rangle\otimes|s_{1}\rangle.$
(1)

Without loss of generality, we may assume that $|\Psi\rangle$ is in canonical
form, so that $M^{s}$ satisfies the condition diagrammatically depicted in Fig
1(a), and by means of the equation

$\sum_{s=1}^{d_{s}}M^{s\dagger}M^{s}=\mathbb{1}_{M},$ (2)

which may be achieved by a suitable definition of $|L\rangle$ and $|R\rangle$
[Pér+07].

Another, more compact, way to represent a MPS (1) is by introducing the
isometry $V:\mathcal{H}_{M}\rightarrow\mathcal{H}_{s}\otimes\mathcal{H}_{M}$,

$V=\sum_{s=1}^{d_{s}}|s\rangle\otimes M^{s},$

satisfying $V^{\dagger}V=\mathbb{1}_{d_{M}}$. With the help of $V$, we may
express (1) as

$|\Psi\rangle=\langle L|\underbrace{V\cdots VV}_{n}|R\rangle,$

where we assume the notational convention that each isometry $V$ acts non-
trivially only on the Hilbert space factor $\mathcal{H}_{M}$, i.e.,
$VV:=(\mathbb{1}_{s}\otimes V)V$. We will henceforth denote the $n$-fold
product $V_{n}:=\underbrace{V\cdots VV}_{n}$. When $|R|$ is used to denote the
number of spins in region $R$, we will write simply $V_{R}$ instead of
$V_{|R|}$.

$\dots$ $\dots$$\dots$ $M$$M$$M$$M$$M$$M$$M$$M$$M$$\dots$ $M$ $M$ $M$ $M$ $M$
$M$ $M$ $M^{\dagger}$$M^{\dagger}$ $M^{\dagger}$ $M^{\dagger}$ $M^{\dagger}$
$M^{\dagger}$ $M^{\dagger}$ $M$$M$ $\tilde{M}_{B}$
$M$$M$$M^{\dagger}$$M^{\dagger}$ $\tilde{M}_{B}^{\dagger}$ $M^{\dagger}$
$M^{\dagger}$$\sigma$$\sigma$$\mathrm{(a)}$$\mathrm{(b)}$$\mathrm{(c)}$$\mathrm{(d)}$$A$$B$$C$
Figure 2: Diagrammatic depictions of: (a) an infinite spin chain with an
example of selection of regions $A$, $B$, and $C$, considered in this paper;
(b) an iuMPS of the spin chain, generated by a tensor $M$ with
$i,j=1,2,...,d_{M}$ and $s=1,2,...,d_{s}$; (c) the reduced density operator
$\rho_{ABC}$ of the iuMPS in (a); the QMC $\tilde{\rho}_{ABC}$ approximating
$\rho_{ABC}$, defined via the tensor $\tilde{M}_{B}$.

An infinite uniform MPS (iuMPS) is a suitable limit of the uniform MPS when
the chain length tends to infinity [FNW92, FNW94].

The limiting state is fully characterized by the collection of reduced density
operators $\rho_{n}$ for a continuous region of $n\in\mathbb{N}$ spins, and
$\rho_{n}$ takes the form

$\displaystyle\rho_{n}$
$\displaystyle=\qquad\sum_{\mathclap{\begin{subarray}{c}s_{1},s_{2},\dots,s_{n}=1\\\
s^{\prime}_{1},s^{\prime}_{2},\dots,s^{\prime}_{n}=1\end{subarray}}}^{d_{s}}\;\mbox{tr}_{M}(M^{s_{n}}\dots
M^{s_{2}}M^{s_{1}}\sigma
M^{s^{\prime}_{1}\dagger}M^{s^{\prime}_{2}\dagger}\dots
M^{s^{\prime}_{n}\dagger})\;|s_{n}\rangle\langle
s^{\prime}_{n}|\otimes\dots\otimes|s_{2}\rangle\langle
s^{\prime}_{2}|\otimes|s_{1}\rangle\langle s^{\prime}_{1}|,$ (3)
$\displaystyle=\mathrm{tr}_{M}(V\cdots VV\sigma V^{\dagger}\cdots
V^{\dagger}V^{\dagger})=\mathrm{tr}_{M}(V_{n}\sigma V^{\dagger}_{n}),$

where $\mbox{tr}\rho_{n}=1,$ and $\sigma$ is a fixed point of the quantum
channel (completely positive trace-preserving (CPTP) map)
$\mathcal{E}:\mathcal{B}(\mathcal{H}_{M})\rightarrow\mathcal{B}(\mathcal{H}_{M})$,
with [FNW92, Pér+07]

$\mathcal{E}(X)=\sum_{s=1}^{d_{s}}M^{s}XM^{s\dagger}.$ (4)

The latter is conveniently represented as a linear operator and defines the
transfer matrix $E\in\mathbb{C}^{d_{M}^{2}\times d_{M}^{2}}$,

$E=\sum_{s=1}^{d_{s}}M^{s}\otimes\bar{M}^{s},$ (5)

where $\bar{M}^{s}$ is the complex conjugate of $M^{s}$. The transfer matrix
is depicted in Fig.1(c), with the matrix elements
$E_{jj^{\prime};ii^{\prime}}:=\langle j|\otimes\langle
j^{\prime}|E|i\rangle\otimes|i^{\prime}\rangle$. To uniquely define a iuMPS we
must specify a family of $d_{M}\times d_{M}$ matrices
$\\{M^{s}\\}_{s=1}^{d_{s}}$ (or, alternatively, the isometry $V$), and the
fixed point density operator $\sigma$ that satisfies
$\mathcal{E}(\sigma)=\sigma$ (refer to Fig. 1(b)).

$\tilde{M}_{B}$$\dots$$i$$j$$s_{1}$$s_{2}$$s_{|B|}$$=:\langle
j|\tilde{M}^{s_{1}s_{2}\dots s_{|B|}}_{B}|i\rangle$ Figure 3: The tensor
$\tilde{M}_{B}$ of the dimension $d_{M}\times d_{s}^{|B|}\times d_{M}$, and
its representation as the matrix elements of the collection of $d_{M}\times
d_{M}$ matrices $\tilde{M}^{s_{1}s_{2}...s_{|B|}}_{B}$, where
$s_{1},s_{2},...,s_{|B|}=1,2,...,d_{s}$.

The eigenvalue problem, depicted in Fig.1(d), and given by

$E|X_{k}\rangle=\nu_{k}|X_{k}\rangle,\quad k=1,\dots,d_{M}^{2},$ (6)

where
$|X_{k}\rangle\in\mathcal{H}_{M}\otimes\mathcal{H}_{M}=\mathbb{C}^{d_{M}}\otimes\mathbb{C}^{d_{M}}$,
determines the spectrum of the transfer matrix. Alternatively, we may
formulate the eigenvalue problem for the quantum channel $\mathcal{E}$,

$\mathcal{E}(X_{k})=\nu_{k}X_{k},\quad k=1,\dots,d_{M}^{2},$

where $X_{k}\in\mathcal{B}(\mathcal{H}_{M})$. For our purposes the most
significant spectral parameter characterizing the transfer matrix is the
spectral gap, characterized by

$\nu_{\mathrm{gap}}=\max_{|\nu_{k}|<1}|\nu_{k}|.$

Note that $\nu_{\mathrm{gap}}<1$ for any $E$, even if several eigenvalues of
$E$ have absolute value one.

Using the notation introduced in (5) we may express $\rho_{n}$ as

$\rho_{n}=\sum_{\mathclap{s_{1},s_{2},\dots,s_{n}=1}}^{d_{s}}\;\langle+|(M^{s_{n}}\otimes\bar{M}^{s^{\prime}_{n}})\cdot\dots\cdot(M^{s_{2}}\otimes\bar{M}^{s^{\prime}_{2}})(M^{s_{1}}\otimes\bar{M}^{s^{\prime}_{1}})|\sigma\rangle\;|s_{n}\rangle\langle
s^{\prime}_{n}|\otimes\dots\otimes|s_{2}\rangle\langle
s^{\prime}_{2}|\otimes|s_{1}\rangle\langle s^{\prime}_{1}|,$ (7)

where $|+\rangle=\sum_{i=1}^{d_{M}}|i\rangle\otimes|i\rangle$ and
$|\sigma\rangle=\sum_{i=1}^{d_{M}}\sigma_{i}|i\rangle\otimes|i\rangle$.

In this paper we divide $n$ spins into three continuous subregions $A$, $B$,
and $C$, and denote the reduced density operator $\rho_{n}$ as $\rho_{ABC}$. A
diagrammatic representation of $\rho_{ABC}$ is given in Fig.2(c). For
$\rho_{ABC}$ we investigate the behavior of $I(A:C|B)$ as a function of the
size of the buffer region $B$, denoted by $|B|$.

##  3 Main Theorem and its proof

Consider the configuration $ABC$ shown in Fig.2(a) and (b), and prepared in an
iuMPS $\rho_{ABC}$, as illustrated in Fig.2(c). For this system the following
theorem on the asymptotic behavior of QCMI holds; it is a tighter version of
the result given in [SK22].

###### Theorem 1

For the described setup of system $ABC$, prepared in an infinite uniform
matrix product state, the quantum conditional mutual information (QCMI) obeys
the bound

$I(A:C|B)\leq Qe^{-q(|B|-K)+2K\ln|B|},$

for $|B|$ large enough. Here $q=-2\ln\nu_{\mathrm{gap}}$, $K$ is a non-
negative integer, and $Q$ does not depend on $|B|$.

Compared to the similar bound in [SK22] with
$q=-\frac{1}{2}\ln\nu_{\mathrm{gap}}$, the value predicted here is
$q=-2\ln\nu_{\mathrm{gap}}.$

Remark: For a given iuMPS $\rho_{ABC}$, the non-negative integer $K$ may be
extracted from the Jordan decomposition of the corresponding transfer matrix
(5) : $K+1$ is the largest dimension of the Jordan blocks corresponding to
eigenvalues of magnitude $\nu_{\mathrm{gap}}$ (see Appendix F for details).

We prove Theorem 1, building on the derivation of Theorem II.1 in [SK22]. For
simplicity, in subsections 3.1 and 3.2 we prove Theorem 1 for the case of
injective iuMPS [Pér+07], i.e., the case for which the transfer matrix (5) has
a unique eigenvector $|\sigma\rangle$ with eigenvalue of magnitude 1. The
generalization to the case of an arbitrary iuMPS is straightforward, but more
technical, and will be discussed in subsection 3.3.

The key idea of the proof is identical to the one used in Theorem II.1 of
[SK22] and involves the construction of a density operator
$\tilde{\rho}_{ABC}$, that is simultaneously a quantum Markov chain and a good
approximation to $\rho_{ABC}$ in trace norm. By comparison with [SK22],
however, we use different methods to estimate $I(A:C|B)$ and
$\|\rho_{ABC}-\tilde{\rho}_{ABC}\|_{1}$, which allow for a tighter bound on
$I(A:C|B)$. The density operator $\tilde{\rho}_{ABC}$ is illustrated
diagrammatically in Fig.2(d) with the related complex tensor $\tilde{M}_{B}$
of dimensions $d_{M}\times d_{s}^{|B|}\times d_{M}$, depicted in Fig.3.

Throughout the Subsections 3.1, 3.2, and 3.3 we exploit the properties of the
the quantum channel
$\tilde{\mathcal{E}}:\mathcal{B}(\mathcal{H}_{M})\rightarrow\mathcal{B}(\mathcal{H}_{M})$
(and its matrix representation
$\tilde{E}:\mathcal{H}_{M}\otimes\mathcal{H}_{M}\rightarrow\mathcal{H}_{M}\otimes\mathcal{H}_{M}$),
which corresponds to the peripheral spectrum of $\mathcal{E}$. It is
straightforward to define $\tilde{\mathcal{E}}$, by first defining $\tilde{E}$
through

$\tilde{E}=PEP,$ (8)

where $P$ is the projector onto the subspace
$\mathrm{span}\\{|u_{i}\rangle\;|\;E|u_{i}\rangle=\nu_{i}|u_{i}\rangle,\,|\nu_{i}|=1\\}$.
The map $\tilde{\mathcal{E}}$ is then defined by

$\langle j|\mathcal{\tilde{E}}(|i\rangle\langle
i^{\prime}|)|j^{\prime}\rangle:=(\langle j|\otimes\langle
j^{\prime}|)\tilde{E}(|i\rangle\otimes|i^{\prime}\rangle).$

It may be shown that the map $\tilde{\mathcal{E}}$ is indeed a quantum channel
(see Lemma III.2 in [SK22]).

The key properties that make $\tilde{E}$ (and $\tilde{\mathcal{E}}$) useful,
are the bounds on $\|E^{n}-\tilde{E}^{n}\|$ derived in Theorem III.2 of
[SRW15] and satisfied for large enough $n$ (specified in Appendix F),

$\displaystyle c_{1}n^{K}\nu_{\mathrm{gap}}^{n-K}$
$\displaystyle\leq\|E^{n}-\tilde{E}^{n}\|$ (9)
$\displaystyle=\|\mathcal{E}^{n}-\tilde{\mathcal{E}}^{n}\|_{2-2}$
$\displaystyle\leq c_{2}n^{K}\nu_{\mathrm{gap}}^{n-K},$

where $n\in\mathbb{N}$, $K$ is a non-negative integer, and $c_{1}$ and $c_{2}$
are some constants, which depend on the details of the Jordan decomposition of
$E$.

###  3.1 A trace distance bound on QCMI.

In this subsection we bound $I(A:C|B)$ in terms of
$\|\rho_{ABC}-\tilde{\rho}_{ABC}\|_{1}$, assuming that $\tilde{\rho}_{ABC}$ is
a quantum Markov chain, which property will be proved in the next subsection.

In QCMI

$\displaystyle I(A:C|B)=\mathrm{tr}$
$\displaystyle\left[\rho_{ABC}\left(\ln\rho_{ABC}+\ln\rho_{B}\right.\right.$
(10) $\displaystyle\left.\left.-\ln\rho_{AB}-\ln\rho_{BC}\right)\right],$

we express

$\rho_{ABC}=\tilde{\rho}_{ABC}+\Delta\rho_{ABC},$

with similar expressions for the marginals of $\rho_{ABC}$ obtained by
suitable partial traces. Then, using the integral representation for the
logarithm of an operator

$\ln\rho=\int_{0}^{\infty}ds\left(\frac{1}{1+s}-\frac{1}{\rho+s}\right),$

we expand (10) to second order in $\Delta\rho_{ABC}$ and its marginals. In
particular for a typical density operator, and suppressing indices, we have

$\displaystyle\ln(\tilde{\rho}+\Delta\rho)=$ $\displaystyle\ln\tilde{\rho}$
(11)
$\displaystyle+\int_{0}^{+\infty}ds\left(\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right.$
$\displaystyle\left.-\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right)$
$\displaystyle+O\left(\|\Delta\rho\|_{1}^{3}\right).$

(Technically, while dealing with logarithm and inverse functions of operators
we should keep track of the subspaces of support. For the sake of clarity, we
proceed as though all density operators $\tilde{\rho}_{ABC}$, etc., are full-
rank, since the rigorous treatment leads to the same final expressions.)

We now consider the zeroth, first and second order terms in $\Delta\rho.$
Since $\tilde{\rho}_{ABC}$ is a quantum Markov chain [LR73a, Rus02],

$\ln\tilde{\rho}_{ABC}+\ln\tilde{\rho}_{B}-\ln\tilde{\rho}_{AB}-\ln\tilde{\rho}_{BC}=0,$

and so the zeroth order term in the expansion of (10) vanishes,

$\displaystyle\mathrm{tr}$
$\displaystyle\left[\tilde{\rho}_{ABC}\left(\ln\tilde{\rho}_{ABC}+\ln\tilde{\rho}_{B}\right.\right.$
$\displaystyle\left.\left.-\ln\tilde{\rho}_{AB}-\ln\tilde{\rho}_{BC}\right)\right]=0.$

The first order contribution has a term that vanishes for the same reason,

$\displaystyle\mathrm{tr}$
$\displaystyle\left[\Delta\rho_{ABC}\left(\ln\tilde{\rho}_{ABC}+\ln\tilde{\rho}_{B}\right.\right.$
$\displaystyle\left.\left.-\ln\tilde{\rho}_{AB}-\ln\tilde{\rho}_{BC}\right)\right]=0.$

The other terms of first order in $\Delta\rho$ are all proportional to an
expression of the following kind,

$\displaystyle\mathrm{tr}\left(\tilde{\rho}\int_{0}^{+\infty}ds\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right)$
$\displaystyle=\mathrm{tr}\left(\int_{0}^{+\infty}ds\frac{1}{\tilde{\rho}+s\mathbb{1}}\tilde{\rho}\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\right)$
$\displaystyle=\mathrm{tr}\left(\left(\int_{0}^{+\infty}ds\frac{1}{(\tilde{\rho}+s\mathbb{1})^{2}}\right)\tilde{\rho}\Delta\rho\right)$
$\displaystyle=\mathrm{tr}\left(\tilde{\rho}^{-1}\tilde{\rho}\Delta\rho\right)$
$\displaystyle=\mathrm{tr}(\Delta\rho)$ $\displaystyle=0,$

since
$\mathrm{tr}\Delta\rho_{ABC}=\mathrm{tr}\rho_{ABC}-\mathrm{tr}\tilde{\rho}_{ABC}=0$
(hence, $\mathrm{tr}\Delta{\rho}_{AB}=0$, etc.).

Thus, the leading non-vanishing terms in the expansion of $I(A:C|B)$ are of
the second order in $\Delta\rho$, coming from the terms of the type

$\displaystyle\mathrm{tr}\left(\Delta\rho\int_{0}^{+\infty}ds\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right)$
(12)
$\displaystyle-\mathrm{tr}\left(\tilde{\rho}\int_{0}^{+\infty}ds\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right)$
$\displaystyle=\mathrm{tr}\left(\int_{0}^{+\infty}ds\frac{s}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right),$

which can be estimated as

$\displaystyle\left|\mathrm{tr}\left(\int_{0}^{+\infty}ds\frac{s}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right)\right|$
(13)
$\displaystyle\leq\int_{0}^{+\infty}dss\left\|\frac{1}{\tilde{\rho}+s\mathbb{1}}\right\|^{3}\|\Delta\rho\|^{2}_{1}$
$\displaystyle\leq\int_{0}^{+\infty}ds\frac{s}{(\lambda_{\mathrm{min}}+s)^{3}}\|\Delta\rho\|^{2}_{1}$
$\displaystyle=\frac{\|\Delta\rho\|^{2}_{1}}{2\lambda_{\mathrm{min}}},$

where $\lambda_{\mathrm{min}}$ is the smallest non-zero eigenvalue of
$\tilde{\rho}$.

Now we estimate $\lambda_{\mathrm{min}}$ with the help of a purification
argument. For the region $R$, where $R$ can be $ABC$, $AB$, $BC$, or $B$, the
density operator $\rho_{R}$ can be purified with the use of the auxiliary
space $\mathcal{H}_{M^{\prime}}$ (note
$\mathcal{H}_{M^{\prime}}\cong\mathcal{H}_{M}$),

$\displaystyle\rho_{R}=\mathrm{tr}_{MM^{\prime}}\left((V_{R}\otimes\mathbb{1}_{M^{\prime}})|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|(V_{R}^{\dagger}\otimes\mathbb{1}_{M^{\prime}})\right),$

where
$|\sqrt{\sigma}\rangle=\sum_{i=1}^{d_{M}}\sqrt{\sigma_{i}}|i\rangle\otimes|i\rangle\in\mathcal{H}_{M}\otimes\mathcal{H}_{M^{\prime}}$,
and $V_{R}$ acts on the space $\mathcal{H}_{M}$. This implies that $\rho_{R}$
has the same spectrum as
$\mathrm{tr}_{R}\left((V_{R}\otimes\mathbb{1}_{M^{\prime}})|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|(V_{R}^{\dagger}\otimes\mathbb{1}_{M^{\prime}})\right)=(\mathcal{E}^{|R|}\otimes\mathrm{id})\left(|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|\right)$,
and hence, assuming
$\dim\mathcal{H}_{R}\geq\dim(\mathcal{H}_{M}\otimes\mathcal{H}_{M^{\prime}})=d_{M}^{2}$,
that

$\rho_{R}=J(\mathcal{E}^{|R|}\otimes\mathrm{id}_{M^{\prime}})\left(|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|\right)J^{\dagger},$
(14)

where
$J:\mathcal{H}_{M}\otimes\mathcal{H}_{M^{\prime}}\rightarrow\mathcal{H}_{R}$
is an isometry. For the purpose of estimating $\lambda_{\mathrm{min}}$ we
define

$\displaystyle\rho^{\prime}_{R}$
$\displaystyle=J(\tilde{\mathcal{E}}^{|R|}\otimes\mathrm{id}_{M^{\prime}})(|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|)J^{\dagger}$
(15)
$\displaystyle=J\left(\sum_{i,j=1}^{d_{M}^{2}}\sqrt{\sigma_{i}\sigma_{j}}\tilde{\mathcal{E}}^{|R|}(|i\rangle\langle
j|)\otimes|i\rangle\langle j|\right)J^{\dagger}.$

The density operator $\rho^{\prime}_{R}$ should not be confused with
$\tilde{\rho}_{R}$. For an injective iuMPS,
$\tilde{\mathcal{E}}:\mathcal{B}(\mathcal{H}_{M})\rightarrow\mathcal{B}(\mathcal{H}_{M})$
acts as $\tilde{\mathcal{E}}(X)=\mathrm{tr}(X)\sigma$, and hence

$\rho^{\prime}_{R}=J\sigma\otimes\sigma J^{\dagger},$ (16)

which implies that the smallest eigenvalue of $\rho^{\prime}_{R}$ is
$\sigma^{2}_{\mathrm{min}},$ where
$\sigma_{\mathrm{min}}:=\mathrm{min}(\sigma_{1},\cdots,\sigma_{d_{M}}).$

To bound $\lambda_{\mathrm{min}}$, we first estimate the trace distance
between $\rho^{\prime}_{R}$ and $\rho_{R}$,

$\displaystyle\|\rho_{R}-\rho^{\prime}_{R}\|_{1}$ (17)
$\displaystyle\quad\leq\|\left((\mathcal{E}^{|R|}-\tilde{\mathcal{E}}^{|R|})\otimes\mathrm{id}\right)(|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|)\|_{1}$
$\displaystyle\quad\leq\|\mathcal{E}^{|R|}-\tilde{\mathcal{E}}^{|R|}\|_{1\rightarrow
1}\||\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|\|_{1}$ $\displaystyle\quad\leq
d_{M}\|\mathcal{E}^{|R|}-\tilde{\mathcal{E}}^{|R|}\|_{2\rightarrow 2}$
$\displaystyle\quad\leq
d_{M}\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow 2}$
$\displaystyle\quad\leq c_{2}d_{M}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$

where we have used the bound (9) in the last line. Then,

$\displaystyle\|\tilde{\rho}_{R}-\rho^{\prime}_{R}\|_{1}$
$\displaystyle\leq\|\tilde{\rho}_{R}-\rho_{R}\|_{1}+\|\rho_{R}-\rho^{\prime}_{R}\|_{1}$
(18) $\displaystyle\leq
c_{2}d_{M}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}+\|\Delta\rho_{R}\|_{1}.$

This implies that,

$\lambda_{\mathrm{min}}=\sigma^{2}_{\mathrm{min}}+O\left(|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\right)+O\left(\|\Delta\rho_{R}\|_{1}\right).$

In the next subsection we will show that $\|\Delta\rho_{R}\|_{1}\leq
O(|B|^{K}\nu_{\mathrm{gap}}^{|B|-K})$. Thus, for $|B|$ large enough,
$\lambda_{\mathrm{min}}$ is (sub)exponentially close to
$\sigma^{2}_{\mathrm{min}}$, and we can guarantee the asymptotic bound

$I(A:C|B)\leq\frac{8}{\sigma^{2}_{\mathrm{min}}}\|\Delta\rho_{ABC}\|_{1}^{2}.$
(19)

This estimate includes the contribution of the remainder term in the expansion
of $I(A:C|B)$ in powers of $\Delta\rho_{ABC}$ that, as shown in Appendix B,
may be bounded by the contribution resulting from terms of the second order.

Remark: The inequality (19) is an upper bound on QCMI quadratic in
$\|\Delta\rho_{ABC}\|_{1}$, in contrast to the upper bound linear in
$\|\Delta\rho_{ABC}\|_{1}$ that we used in [SK22]. In Subsection 3.2 we will
see that this improved bound is responsible for changing the rate of
exponential decay of the QCMI bounding function by a factor of $2$, compared
to that in Theorem II.1 in [SK22].

###  3.2 A (sub)exponential bound on trace distance.

In this subsection we implicitly construct the tensor $\tilde{M}_{B}$, such
that $\tilde{\rho}_{ABC}$ defined in Fig.2 becomes a quantum Markov chain, and
we show that

$\displaystyle\|\Delta\rho_{ABC}\|_{1}$ $\displaystyle\leq
O(\|E^{|B|}-\tilde{E}^{|B|}\|)$ (20) $\displaystyle\leq
O(|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}),$

as a result. We represent the tensor $\tilde{M}_{B}$ as a collection of
$d_{s}^{|B|}$ complex $d_{M}\times d_{M}$ matrices $\tilde{M}^{s_{1}s_{2}\dots
s_{|B|}}_{B}$, where $s_{k}=1,2,\dots,d_{s}$ for $k=1,2,\dots,|B|$. The
matrices $\tilde{M}^{s_{1}s_{2}\dots s_{|B|}}_{B}$ are defined by their matrix
elements in Fig. 3.

For an injective iuMPS, the quantum channel $\tilde{\mathcal{E}}$ has the form
[FNW92]

$\tilde{\mathcal{E}}(X)=\mathrm{tr}(X)\sigma.$

Its matrix representation $\tilde{E}$ is then obtained with the use of (5),

$\tilde{E}=\tilde{E}^{n}=|\sigma\rangle\langle+|.$

In Section IV.A of [SK22] it was proved, that if we impose on $\tilde{M}_{B}$
the condition

$\displaystyle\quad\;\mathclap{\sum_{s_{1},\dots,s_{|B|}=1}^{d_{s}}}\quad\tilde{M}_{B}^{s_{1}\dots
s_{|B|}}\otimes\overline{\tilde{M}}_{B}^{s_{1}\dots
s_{|B|}}=\tilde{E}^{|B|}=|\sigma\rangle\langle+|,$ (21)

or, equivalently, the condition

$\displaystyle\tilde{\mathcal{E}}^{|B|}(X)$
$\displaystyle=\sum_{s_{1},\dots,s_{|B|}=1}^{d_{s}}\tilde{M}_{B}^{s_{1}\dots
s_{|B|}}X\tilde{M}_{B}^{s_{1}\dots s_{|B|}\dagger}$ (22)
$\displaystyle=\mathrm{tr}(X)\sigma,$

then the induced $\tilde{\rho}_{ABC}$ is a quantum Markov chain. The
$d_{M}\times d_{M}$ complex matrices $\overline{\tilde{M}}_{B}^{s_{1}\dots
s_{|B|}}$ and $\tilde{M}_{B}^{s_{1}\dots s_{|B|}\dagger}$ are the element-wise
complex conjugate and the complex conjugate transpose of
$\tilde{M}_{B}^{s_{1}\dots s_{|B|}}$, respectively. We will now show that the
conditions (21) and (22) result in the desired bound (20).

Notice that the isometries

$\displaystyle V_{B}$
$\displaystyle=\sum_{s_{1},\dots,s_{|B|}=1}^{d_{s}}|s_{1}\rangle\otimes|s_{2}\rangle\otimes\cdots\otimes|s_{|B|}\rangle\otimes
M^{s_{1}}M^{s_{2}}\dots M^{s_{|B|}},$ $\displaystyle\tilde{V}_{B}$
$\displaystyle=\sum_{s_{1},\dots,s_{|B|}=1}^{d_{s}}|s_{1}\rangle\otimes|s_{2}\rangle\otimes\cdots\otimes|s_{|B|}\rangle\otimes\tilde{M}_{B}^{s_{1}\dots
s_{|B|}},$

represent Stinespring dilations [Sti55] for the quantum channels
$\mathcal{E}^{|B|}$ and $\tilde{\mathcal{E}}^{|B|}$, respectively, with
$\mathcal{H}_{B}$ serving as the dilating Hilbert space, i.e.,

$\displaystyle\mathcal{E}^{|B|}(X)$
$\displaystyle=\mathrm{tr}_{B}(V_{B}XV_{B}^{\dagger}),$ (23)
$\displaystyle\tilde{\mathcal{E}}^{|B|}(X)$
$\displaystyle=\mathrm{tr}_{B}(\tilde{V}_{B}X\tilde{V}_{B}^{\dagger})=\mathrm{tr}(X)\sigma.$
(24)

Moreover, we can conveniently reexpress $\tilde{M}_{B}$ in terms of the
isometry $\tilde{V}_{B}$ as

$\tilde{M}_{B}^{s_{1}\dots s_{|B|}}=(\langle s_{1}|\otimes\langle
s_{2}|\otimes\cdots\otimes\langle s_{|B|}|)\tilde{V}_{B},$ (25)

which allows us to work with $\tilde{V}_{B}$ satisfying (23) instead of
$\tilde{M}_{B}$ satisfying (22).

To bound $\|\Delta\rho_{ABC}\|_{1}$, we notice that we can express
$\rho_{ABC}$ and $\tilde{\rho}_{ABC}$ in terms of $V_{B}$ and $\tilde{V}_{B}$
as

$\displaystyle\rho_{ABC}$ $\displaystyle=\mathrm{tr}_{M}(V_{C}V_{B}V_{A}\sigma
V^{\dagger}_{A}V_{B}^{\dagger}V^{\dagger}_{C}),$
$\displaystyle\tilde{\rho}_{ABC}$
$\displaystyle=\mathrm{tr}_{M}(V_{C}\tilde{V}_{B}V_{A}\sigma
V^{\dagger}_{A}\tilde{V}_{B}^{\dagger}V^{\dagger}_{C}).$

It follows that [SK22],

$\displaystyle\|\Delta\rho_{ABC}\|_{1}\leq 2\|V_{B}-\tilde{V}_{B}\|.$ (26)

Since the dimension of the dilating Hilbert space
$\dim\mathcal{H}_{B}=d_{s}^{|B|}$ is exponentially large, an explicit
construction of $\tilde{V}_{B}$ is not useful. Instead we notice that the
quantum channels $\mathcal{E}^{|B|}$ and $\tilde{\mathcal{E}}^{|B|}$ may be
dilated with a Hilbert space of dimension $d_{M}^{2}$, i.e. a space isomorphic
to $\mathcal{H}_{M}\otimes\mathcal{H}_{M}$. We will thus explicitly construct
the isometries
$W,\tilde{W}:\mathcal{H}_{M}\rightarrow\mathcal{H}_{M}\otimes\mathcal{H}_{M}\otimes\mathcal{H}_{M}$,
representing Stinespring dilations of $\mathcal{E}^{|B|}$ and
$\tilde{\mathcal{E}}^{|B|}$, respectively. We may then express $V_{B}$ and
$\tilde{V}_{B}$ in terms of these isometries, as we now show. Since $W$ and
$V_{B}$ represent dilations of the same quantum channel $\mathcal{E}^{|B|}$,
then, by Stinespring’s theorem [Sti55] and assuming $\dim\mathcal{H}_{B}\geq
d_{M}^{2}$, _there exists_ an isometry
$U:\mathcal{H}_{M}\otimes\mathcal{H}_{M}\rightarrow\mathcal{H}_{B}$, such that

$V_{B}=(U\otimes\mathbb{1}_{M})W.$

Using this isometry $U$, we may re-express $\tilde{V}_{B}$ as

$\tilde{V}_{B}:=(U\otimes\mathbb{1}_{M})\tilde{W}.$ (27)

We will be able to construct $\tilde{W}$ explicitly, and thus $\tilde{V}_{B}$
follows from (27) and $\tilde{M}_{B}$ follows from (25).

We may now express the bound (26) in terms of the norm distance between $W$
and $\tilde{W}$,

$\displaystyle\|\Delta\rho_{ABC}\|_{1}$ $\displaystyle\leq
2\|V_{B}-\tilde{V}_{B}\|$ (28)
$\displaystyle=2\|(U\otimes\mathbb{1}_{M})(W-\tilde{W})\|$
$\displaystyle=2\|W-\tilde{W}\|,$

eliminating the need to explicitly specify $U$, which we will not do.

We now construct $W$ and $\tilde{W}$ using Choi’s methods [Cho75], which will
enable us to show that $\|\Delta\rho_{ABC}\|_{1}\leq
O(\|E^{|B|}-\tilde{E}^{|B|}\|)$.

Let $C$ and $\tilde{C}$ be the Choi matrices of $\mathcal{E}^{|B|}$ and
$\tilde{\mathcal{E}}^{|B|}$, respectively,

$\displaystyle C$
$\displaystyle=\left(\mathrm{id}\otimes\mathcal{E}^{|B|}\right)\left(|+\rangle\langle+|\right),$
(29) $\displaystyle\tilde{C}$
$\displaystyle=\left(\mathrm{id}\otimes\tilde{\mathcal{E}}^{|B|}\right)\left(|+\rangle\langle+|\right)=\sigma\otimes\mathbb{1}_{M},$

where we have used (23). We further define

$\Delta
C:=C-\tilde{C}=\left(\mathrm{id}\otimes\left(\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\right)\right)\left(|+\rangle\langle+|\right),$
(30)

from which it follows

$\displaystyle\|\Delta C\|_{1}$ $\displaystyle\leq
d_{M}^{2}\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow 2}$
(31) $\displaystyle\leq c_{2}d_{M}^{2}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$
$\displaystyle\|\Delta C\|$ $\displaystyle\leq
d_{M}\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow 2}$
$\displaystyle\leq c_{2}d_{M}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$

where we have used (9). We aim to express $\|W-\tilde{W}\|$ in terms of
$\|\Delta C\|_{1}$ and $\|\Delta C\|$, obtaining the desired bound (20).

We now employ Choi’s construction [Cho75]. Recall that
$|+\rangle:=\sum_{i=1}^{d_{M}}|i\rangle\otimes|i\rangle$ in (29) is defined in
our chosen orthonormal product basis, which satisfies
$(\mathbb{1}\otimes\sigma)|i\rangle\otimes|j\rangle=\sigma_{j}|i\rangle\otimes|j\rangle$
and
$(\sigma\otimes\mathbb{1})|i\rangle\otimes|j\rangle=\sigma_{i}|i\rangle\otimes|j\rangle$
for all $i,j=1,2,\dots,d_{M}$. With respect to these basis, we decompose

$\displaystyle\tilde{C}=\sum_{i,j=1}^{d_{M}}|\tilde{\psi}_{ij}\rangle\langle\tilde{\psi}_{ij}|,$
(32)
$\displaystyle|\tilde{\psi}_{ij}\rangle:=\sqrt{\sigma_{i}}|i\rangle\otimes|j\rangle.$

Notice that the vectors $|\tilde{\psi}_{ij}\rangle$ are pairwise orthogonal,
$\langle\tilde{\psi}_{i^{\prime}j^{\prime}}|\tilde{\psi}_{ij}\rangle=\sqrt{\sigma_{i}\sigma_{i^{\prime}}}\delta_{ii^{\prime}}\delta_{jj^{\prime}}$,
and form a basis in $\mathcal{H}_{M}\otimes\mathcal{H}_{M}$. We define a set
of projections
$P_{k}:\mathcal{H}_{M}\otimes\mathcal{H}_{M}\rightarrow\mathcal{H}_{M}$ with
$k=1,2,\dots,d_{M}$ by

$P_{k}:=\mathbb{1}_{M}\otimes\langle k|.$ (33)

Now Choi’s construction [Cho75] provides the isometry $\tilde{W}$ dilating
$\tilde{\mathcal{E}}^{|B|}$ explicitly in the form

$\displaystyle\tilde{W}:=\sum_{k,i,j=1}^{d_{M}}\left(P_{k}|\tilde{\psi}_{ij}\rangle\otimes|i\rangle\otimes|j\rangle\right)\langle
k|.$

Now we assume a decomposition of $C$ of the form,

$C=\sum_{i,j=1}^{d_{M}}|\psi_{ij}\rangle\langle\psi_{ij}|,$ (34)

where the vectors $|\psi_{ij}\rangle$ are not necessarily orthogonal. Recall
that Choi’s construction does not require the vectors $|\psi_{ij}\rangle$ in
the decomposition (34) to be orthogonal. Using the projectors $P_{k}$ defined
in (33), we employ Choi’s construction again, obtaining the isometry

$\displaystyle
W:=\sum_{k,i,j=1}^{d_{M}}\left(P_{k}|\psi_{ij}\rangle\otimes|i\rangle\otimes|j\rangle\right)\langle
k|,$ (35)

which dilates $\mathcal{E}^{|B|}$.

Since $|\tilde{\psi}_{ij}\rangle$ form a basis, we can express,

$\displaystyle|\psi_{ij}\rangle$
$\displaystyle=|\tilde{\psi}_{ij}\rangle+\sum_{m,n=1}^{d_{M}}\gamma_{mn;ij}|\tilde{\psi}_{mn}\rangle$
(36) $\displaystyle:=\left(\mathbb{1}_{M\otimes
M}+\Gamma\right)|\tilde{\psi}_{ij}\rangle,$

where $\gamma_{mn;ij}\in\mathbb{C}$ are some coefficients, and where we have
defined

$\Gamma:=\sum_{i,j,m,n=1}^{d_{M}}\gamma_{mn;ij}\sigma^{-1}_{i}|\tilde{\psi}_{mn}\rangle\langle\tilde{\psi}_{ij}|.$

We express the distance in operator norm between the constructed isometries
$W$ and $\tilde{W}$ in terms of $\Gamma$,

$\displaystyle\|W-\tilde{W}\|^{2}$
$\displaystyle=\|(W-\tilde{W})^{\dagger}(W-\tilde{W})\|$
$\displaystyle\leq\|(W-\tilde{W})^{\dagger}(W-\tilde{W})\|_{1}$
$\displaystyle=\mathrm{tr}\left((W-\tilde{W})^{\dagger}(W-\tilde{W})\right)$
$\displaystyle=\sum_{i,j,m,n=1}^{d_{M}}|\gamma_{mn;ij}|^{2}\sigma_{m}$
$\displaystyle=\mathrm{tr}\left(\Gamma(\sigma\otimes\mathbb{1}_{M})\Gamma^{\dagger}\right)$
$\displaystyle=\mathrm{tr}\left(\Gamma\tilde{C}\Gamma^{\dagger}\right).$

We will construct such $\Gamma$, that (34) is consistent with (36), and so
that $\|W-\tilde{W}\|=O(\|\Gamma\|_{1})=O(\|\Delta C\|_{1})$. To find
$\Gamma$, we proceed as follows, and note that from (32), (34), and (36),

$\displaystyle\Delta C$ $\displaystyle=C-\tilde{C}$ (37)
$\displaystyle=\left(\mathbb{1}_{M\otimes
M}+\Gamma\right)\tilde{C}(\mathbb{1}_{M\otimes M}+\Gamma^{\dagger})-\tilde{C}$
$\displaystyle=\Gamma\tilde{C}+\tilde{C}\Gamma^{\dagger}+\Gamma\tilde{C}\Gamma^{\dagger}.$

We make the definitions,

$\displaystyle\Delta C^{\prime}$
$\displaystyle:=\tilde{C}^{-\frac{1}{2}}\Delta C\tilde{C}^{-\frac{1}{2}}$ (38)
$\displaystyle=(\sigma^{-\frac{1}{2}}\otimes\mathbb{1}_{M})\Delta
C(\sigma^{-\frac{1}{2}}\otimes\mathbb{1}_{M}),$ $\displaystyle\Gamma^{\prime}$
$\displaystyle:=\tilde{C}^{-\frac{1}{2}}\Gamma\tilde{C}^{\frac{1}{2}}$
$\displaystyle=(\sigma^{-\frac{1}{2}}\otimes\mathbb{1}_{M})\Gamma(\sigma^{\frac{1}{2}}\otimes\mathbb{1}_{M}),$

transforming the equation (37) into

$\Gamma^{\prime}+\Gamma^{\prime\dagger}+\Gamma^{\prime}\Gamma^{\prime\dagger}=\Delta
C^{\prime}.$ (39)

Notice that from (38) it follows that

$\|\Delta C^{\prime}\|_{1}\leq\frac{1}{\sigma_{\mathrm{min}}}\|\Delta
C\|_{1}\leq\frac{c_{2}d^{2}_{M}}{\sigma_{\mathrm{min}}}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$
(40)

where we have used (31). This implies that $\Delta C^{\prime}\rightarrow 0$ as
$|B|\rightarrow\infty$. We wish to find $\Gamma^{\prime}$ satisfying (39),
such that $\|\Gamma^{\prime}\|_{1}=O(\|\Delta C^{\prime}\|_{1})$, and via
(40), $\|\Gamma\|_{1}=O(\|\Delta C\|_{1})$. Since any solution of (39)
satisfying this condition will suffice, and observing that $\Delta C^{\prime}$
is Hermitian, we shall seek a Hermitian solution,
$\Gamma^{\prime}=\Gamma^{\prime\dagger}$.

For $\Gamma^{\prime}$ Hermitian, equation (39) reduces to a quadratic
equation,

$\Gamma^{\prime 2}+2\Gamma^{\prime}-\Delta C^{\prime}=0,$

which indeed has a Hermitian solution of small norm,

$\Gamma^{\prime}=\sqrt{\mathbb{1}_{M\otimes M}+\Delta
C^{\prime}}-\mathbb{1}_{M\otimes M},$ (41)

provided that $\mathbb{1}_{M\otimes M}+\Delta C^{\prime}\geq 0$, which is
satisfied for large enough $|B|$. Notice that through (38), (36), and (35),
$\Gamma^{\prime}$ explicitly defines $W$.

From (41) follows,

$\Gamma^{\prime}\Gamma^{\prime\dagger}=\Gamma^{\prime 2}=2\mathbb{1}_{M\otimes
M}+\Delta C^{\prime}-2\sqrt{\mathbb{1}_{M\otimes M}+\Delta C^{\prime}}.$

Assuming $|B|$ large enough that $\|\Delta C^{\prime}\|\leq 2^{-\frac{3}{2}}$,
by Taylor expansion, we estimate

$0\leq\Gamma^{\prime 2}\leq\frac{\Delta C^{\prime 2}}{2}.$ (42)

Expressing $\|W-\tilde{W}\|$ in terms of $\Gamma^{\prime}$, we find

$\|W-\tilde{W}\|^{2}=\mathrm{tr}\left(\Gamma\tilde{C}\Gamma^{\dagger}\right)=\mathrm{tr}\left(\tilde{C}\Gamma^{\prime}\Gamma^{\prime\dagger}\right),$

leading to

$\displaystyle\|W-\tilde{W}\|^{2}$
$\displaystyle\leq\frac{1}{2}\mathrm{tr}\left(\tilde{C}\Delta C^{\prime
2}\right)$ (43) $\displaystyle=\frac{1}{2}\mathrm{tr}\left(\Delta
C^{2}\tilde{C}^{-1}\right)$ $\displaystyle\leq\frac{1}{2\sigma_{\min}}\|\Delta
C\|\|\Delta C\|_{1}$
$\displaystyle\leq\frac{d_{M}^{3}}{2\sigma_{\min}}\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow
2}^{2},$

where we have used (29), (31), and (38).

Combining (26), (28), (43), and (9) together, we have

$\|\Delta\rho_{ABC}\|_{1}\leq
c_{2}d_{M}\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$
(44)

proving that (an injective) iuMPS is an approximate quantum Markov chain.

Finally, using (19), we obtain the desired bound on QCMI,

$I(A:C|B)\leq\frac{16d_{M}^{3}c_{2}^{2}}{\sigma_{\min}^{3}}e^{-2\ln\nu^{-1}_{\mathrm{gap}}\,(|B|-K)+2K\ln|B|},$
(45)

which corresponds to the statement of Theorem 1 with
$Q=\frac{16d_{M}^{3}c_{2}^{2}}{\sigma_{\min}^{3}}$,
$q=2\ln\nu^{-1}_{\mathrm{gap}}$. This proves Theorem (1) for the case of
injective iuMPS.

Remarks:

  1. 1.

Notice that the inequality (44) has the form $\|\Delta\rho_{ABC}\|_{1}\leq
O(\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow 2})$, in
contrast to the bound $\|\Delta\rho_{ABC}\|_{1}\leq
O(\sqrt{\|\mathcal{E}^{|B|}-\tilde{\mathcal{E}}^{|B|}\|_{2\rightarrow 2}})$
that we have obtained in Theorem II.1 of [SK22]. This improvement modifies the
rate of exponential decay in (45) by a factor of $2$, in addition to the
factor of $2$ discussed in Remark 1 in Subsection 3.1. The rate of the
exponential decay in the final bound (45) is thus four times larger than the
one in Theorem II.1 of [SK22].

  2. 2.

The numerical results presented in Section 4 are much closer to the new bound
on the decay rate than the earlier value given in [SK22].

  3. 3.

Another improvement over Theorem II.1 of [SK22] is that the prefactor in (45)
is independent of the size of any region of the spin chain. As a drawback, we
gained the dependence on $\sigma^{-3}_{\mathrm{min}}$, which may be a large
number.

###  3.3 General case: infinite uniform matrix product states.

Now we outline the changes required to extend the proof to the general class
of iuMPS. While the main steps of the proof remain the same, we have to
account for the cases in which iuMPS is not injective. These entail changes in
the peripheral part $\tilde{E}$ (8) of the transfer matrix $E$ (5). Since the
modifications we make involve the equations (15) and (29), formulated in terms
of the quantum channel $\tilde{\mathcal{E}}$ (represented by $\tilde{E}$), we
focus attention on $\tilde{\mathcal{E}}$, rather than $\tilde{E}$.

In the case of a general iuMPS, the quantum channel
$\tilde{\mathcal{E}}:\mathcal{B}(\mathcal{H}_{M})\rightarrow\mathcal{B}(\mathcal{H})_{M}$
has the following properties [FNW92] (see also Propositions 1 and 2 in
[SK22]):

  * •

The quantum channel $\tilde{\mathcal{E}}$ decomposes into a sum of components,
$\tilde{\mathcal{E}}=\sum_{\alpha=1}^{N}\tilde{\mathcal{E}}_{\alpha}$,
$\tilde{\mathcal{E}}_{\alpha}:\Pi^{\alpha}\mathcal{B}(\mathcal{H}_{M})\Pi^{\alpha}\rightarrow\Pi^{\alpha}\mathcal{B}(\mathcal{H}_{M})\Pi^{\alpha}$,
defined by the set of projectors $\\{\Pi^{\alpha}\\}_{\alpha=1}^{N}$,
$\sum_{\alpha=1}^{N}\tilde{\Pi}^{\alpha}=\mathbb{1}_{M}$. This implies that
the fixed point $\sigma$ of $\mathcal{E}$ and $\tilde{\mathcal{E}}$ satisfies
$\sigma=\sum_{\alpha=1}^{N}\Pi^{\alpha}\sigma\Pi^{\alpha}$.

  * •

For each $\alpha$, the component $\tilde{\mathcal{E}}_{\alpha}$ admits a
further decomposition,

$\tilde{\mathcal{E}}_{\alpha}(X)=p_{\alpha}\sum_{k=0}^{p_{\alpha}-1}\mathrm{tr}\left(\Pi^{\alpha}_{k}X\Pi^{\alpha}_{k}\right)(\Pi^{\alpha}_{k+1}\sigma\Pi^{\alpha}_{k+1}),$

defined by the set of projectors $\Pi^{\alpha}_{k}$ with
$k=0,\dots,p_{\alpha}-1$, satisfying
$\sum_{k=0}^{p_{\alpha}-1}\Pi^{\alpha}_{k}=\Pi^{\alpha}$. The index $k$ in
$\Pi^{\alpha}_{k}$ is taken modulo $p_{\alpha}$. This implies that $\sigma$
satisfies
$\sigma=\sum_{\alpha=1}^{N}\sum_{k=0}^{p_{\alpha}-1}\Pi_{k}^{\alpha}\sigma\Pi_{k}^{\alpha}$.

We first consider the case when $N=1$, where there is a single component in
the direct sum decomposition
$\tilde{\mathcal{E}}=\bigoplus_{\alpha=1}^{N}\tilde{\mathcal{E}}_{\alpha}$. We
allow, however, arbitrary period $p\geq 1$. Under these assumptions and using
(29),

$\tilde{C}=p\sum_{k=0}^{p-1}\Pi_{k+|B|}\sigma\Pi_{k+|B|}\otimes\Pi_{k},$

which implies

$\|\tilde{C}^{-1}\|=\frac{1}{p\sigma_{\min}}.$

This modifies the pre-factor in the bound (43) from
$\frac{d_{M}^{3}}{2\sigma_{\min}}$ to $\frac{d_{M}^{3}}{2p\sigma_{\min}}$.

The density operator $\rho^{\prime}_{R}$, defined in (15), should be updated
from (16) to

$\displaystyle\rho^{\prime}_{R}=pJ\sum_{k=0}^{p-1}\Pi_{k}\sigma\Pi_{k}\otimes\sqrt{\sigma}\Pi_{k+|R|}\sqrt{\sigma}J^{\dagger}$
$\displaystyle=pJ\sum_{k=0}^{p-1}\Pi_{k}\sigma\Pi_{k}\otimes\Pi_{k+|R|}\sigma\Pi_{k+|R|}J^{\dagger}$
$\displaystyle=pJ\sum_{k=0}^{p-1}\left(\Pi_{k}\otimes\Pi_{k+|R|}\right)(\sigma\otimes\sigma)$
$\displaystyle\qquad\qquad\qquad\qquad\times\left(\Pi_{k}\otimes\Pi_{k+|R|}\right)J^{\dagger},$

where $J:\mathcal{H}_{M}\otimes\mathcal{H}_{M}\rightarrow\mathcal{H}_{R}$ is
an isometry, and where we have used
$\sigma=\sum_{k=0}^{p-1}\Pi_{k}\sigma\Pi_{k}$. This implies that for any
$\rho_{R}$ or $\tilde{\rho}_{R}$, via (17) and (18), its smallest eigenvalue
$\lambda_{\mathrm{min}}$ is bounded from below with exponential precision by

$\lambda_{\mathrm{min}}\geq
p\sigma^{2}_{\mathrm{min}}+O\left(|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\right).$

This modifies the pre-factor in the bound (19) from
$\frac{8}{\sigma^{2}_{\mathrm{min}}}$ to
$\frac{8}{p\sigma^{2}_{\mathrm{min}}}$. Thus, the overall bound (45) gains a
factor $p^{-2}\leq 1$,

$I(A:C|B)\leq\frac{16d_{M}^{3}c_{2}^{2}}{p^{2}\sigma_{\min}^{3}}e^{-2\ln\nu^{-1}_{\mathrm{gap}}\,(|B|-K)+2K\ln|B|},$
(46)

proving Theorem 1 for the case $N=1$ and $p\geq 1$.

Generalizing to the case of $N\geq 1$ we notice that $\tilde{C}$ obtains a
direct sum structure,

$\tilde{C}=\bigoplus_{\alpha=1}^{N}p_{\alpha}\sum_{k=0}^{p_{\alpha}-1}\Pi^{\alpha}_{k+|B|}\sigma\Pi^{\alpha}_{k+|B|}\otimes\Pi^{\alpha}_{k}.$
(47)

The density operator $\rho^{\prime}_{R}$ also gets a direct sum structure,

$\rho^{\prime}_{R}=J\bigoplus_{\alpha=1}^{N}p_{\alpha}\sum_{k=0}^{p_{\alpha}-1}\left(\Pi^{\alpha}_{k}\otimes\Pi^{\alpha}_{k+|R|}\right)\\\
(\sigma\otimes\sigma)\left(\Pi^{\alpha}_{k}\otimes\Pi^{\alpha}_{k+|R|}\right)J^{\dagger}.$
(48)

Proceeding as in the case of $N=1$, the equations (47) and (48) imply that the
bound (44) is modified to

$\displaystyle\|\Delta\rho_{ABC}\|_{1}\leq$ (49)
$\displaystyle\qquad\max_{\alpha}\left(\sqrt{\frac{d^{3}_{M_{\alpha}}}{p_{\alpha}}}\right)c_{2}\sqrt{\frac{2}{\sigma_{\mathrm{min}}}}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K},$

and the bound (46) to

$\displaystyle I(A:C|B)\leq$ (50)
$\displaystyle\quad\max_{\alpha}\left(\frac{d_{M_{\alpha}}^{3}}{p_{\alpha}^{2}}\right)\frac{16c_{2}^{2}}{\sigma_{\min}^{3}}e^{-2\ln\nu^{-1}_{\mathrm{gap}}\,(|B|-K)+2K\ln|B|},$

where $d_{M_{\alpha}}=\dim(\Pi^{\alpha}\mathcal{H}_{M})$. The bound (49)
implies that any iuMPS is an approximate quantum Markov chain for sufficiently
large $|B|$. We note that since $p_{\alpha}\geq 1$ and $d_{M_{\alpha}}\leq
d_{M}$, the bounds (44) and (45) hold for the general class of iuMPS.

##  4 Numerical calculation of QCMI

###  4.1 Constructing iuMPS

$M$ $M$ $\dots$ $M$ $M$ $M^{\dagger}$ $M^{\dagger}$ $\dots$ $M^{\dagger}$
$M^{\dagger}$
$i$$i^{\prime}$$j$$j^{\prime}$$=$$G$$=$$\Sigma$$W$$W^{\dagger}$$j$$i$$i^{\prime}$$j^{\prime}$$l$$i$$j$$l^{\prime}$$i^{\prime}$$j^{\prime}$
Figure 4: The first equality in the diagram depicts the definition of
$G_{ii^{\prime};jj^{\prime}}$; the second equality depicts its
eigendecomposition.

Recall that an iuMPS is completely characterized by reduced density operators
$\rho_{n}$ on finite continuous regions of $n$ sites (refer to (7)). In turn,
$\rho_{n}$ is characterized by the fixed point of the quantum channel
$\mathcal{E}$ and the choice of the set $\mathcal{S}$. The quantum channel
$\mathcal{E}$ is represented by a transfer matrix $E$ that may have a single
or multiple fixed point(s) (eigenvectors corresponding to eigenvalue 1 of
$E$). There are three cases that account for all iuMPS:

Case 1. When $E$ has a single fixed point, and,

$\displaystyle\mathcal{S}:=\left\\{\left.M^{s}\right\rvert
M^{s}\in\mathbb{C}^{d_{M}\times d_{M}},\forall
s\in\left\\{1,2,\dots,d_{s}\right\\}\right\\},$

Case 2. When $E$ has two fixed points and $\rho_{n}$ is of period 1, and,

$\mathcal{S}:=\\\ \left\\{\left.M^{s}:=\begin{bmatrix}M_{1}^{s}&0\\\
0&M_{2}^{s}\end{bmatrix}\right\rvert
M_{1}^{s},M_{2}^{s}\in\mathbb{C}^{\frac{d_{M}}{2}\times\frac{d_{M}}{2}},\right.\\\
\left.M^{s}\in\mathbb{C}^{d_{M}\times d_{M}},\forall
s\in\left\\{1,2,\dots,d_{s}\right\\}\right.\Bigg{\\}},$

Case 3. When $E$ has two fixed points and $\rho_{n}$ is of period 2,

$\mathcal{S}:=\\\ \left\\{\left.M^{s}:=\begin{bmatrix}0&M_{1}^{s}\\\
M_{2}^{s}&0\end{bmatrix}\right\rvert
M_{1}^{s},M_{2}^{s}\in\mathbb{C}^{\frac{d_{M}}{2}\times\frac{d_{M}}{2}},\right.\\\
\left.M^{s}\in\mathbb{C}^{d_{M}\times d_{M}},\forall
s\in\left\\{1,2,\dots,d_{s}\right\\}\right.\Bigg{\\}},$

The matrices $M^{s}$ in the set $\mathcal{S}$ in Case 1 are constructed as
follows,

$\displaystyle\left(M^{s}\right)_{ij}=\sum_{s^{\prime}=1}^{d_{s}}\sum_{j^{\prime}=1}^{d_{M}}U_{ij^{\prime}}^{ss^{\prime}}\Psi_{j^{\prime}j}^{s^{\prime}},$
(51)

where the unitary $U\in\mathbb{C}^{d_{s}d_{M}\times d_{s}d_{M}}$ is selected
randomly with respect to Haar measure and $\Psi\in\mathbb{C}^{d_{s}d_{M}\times
d_{M}}$ is given by,

$\displaystyle\Psi:=\frac{1}{\sqrt{3}}\begin{bmatrix}1\\\ 1\\\
1\end{bmatrix}\otimes\mathds{1}_{d_{M}},$ (52)

where we fix $d_{s}=3$ and $d_{M}=4$ in our numerical calculations. In Case 2
and Case 3, each of the two matrices $M_{1}^{s}$ and $M_{2}^{s}$ is
constructed using an equation of the form (51) and (52), but with $d_{M}$
replaced by $\frac{d_{M}}{2}$, and different and independent choices of the
unitary $U$.

After constructing $\mathcal{S}$, the left and right fixed points of the
transfer matrix $E$ denoted by $|+\rangle$ and $|\sigma\rangle$ in (7),
respectively, remain to be found to completely specify the iuMPS. Recall,
fixed points of $E$ are its left and right eigenvectors corresponding to a
left and right eigenvalue of 1, respectively. If there are multiple fixed
points of $E$, in other words, if (left and/or right) eigenvalue of 1 of $E$
is degenerate, then the fixed point used to define iuMPS is the uniform linear
combination of the multiple fixed points. That is, the left and/or right
eigenvectors are defined as the normalized sum of all left and/or right
eigenvectors, respectively.

###  4.2 Calculating entropies

$\Sigma^{-\frac{1}{2}}$$W^{\dagger}$$M$ $M$ $\dots$ $M$ $M$ $M^{\dagger}$
$M^{\dagger}$ $\dots$ $M^{\dagger}$ $M^{\dagger}$
$W$$\Sigma^{-\frac{1}{2}}$$=$$\Sigma$$W$$W^{\dagger}$$\Sigma^{-\frac{1}{2}}$$W^{\dagger}$$W$$\Sigma^{-\frac{1}{2}}$$l$$l^{\prime}$$l$$l^{\prime}$$=$$l$$l^{\prime}$
Figure 5: Verification of the property $P^{\dagger}P=\mathds{1}_{d_{M}^{2}}$.
$\sigma$$W^{\dagger}$$\Sigma^{\frac{1}{2}}$$W$$\Sigma^{\frac{1}{2}}$$l$$l^{\prime}$
Figure 6: The diagram depicts $P^{\dagger}\rho_{n}P$, the density operator
$\rho_{n}$ projected onto its support, where $P$ is defined by (55) and
satisfies the property
$P^{\dagger}\mathds{1}_{d_{s}^{n}}P=\mathds{1}_{d_{M}^{2}}$.

QCMI is computed from von Neumann entropies of the states on the continuous
subregions $AB$, $BC$, $B$ and $ABC$ that are in turn computed from the
eigenvalues of the reduced density operators for the respective subregions. A
numerical check of Theorem 1 requires QCMI for increasing size, $\left\lvert
B\right\rvert$, of the separating region $B$. Direct diagonalization is
computationally inefficient, because even linear growth of $\left\lvert
B\right\rvert$ would imply exponential growth of the Hilbert space dimension
of any subregion that includes $B$. To overcome this fundamental issue of
dimensionality, we use tensor network contractions to project reduced density
operators on to their respective supports222The benefit of finding eigenvalues
of $\rho_{n}$ by projecting $\rho_{n}$ onto its support versus direct
diagonalization is revealed by comparing how the memory required to implement
the two schemes scales with $n$. Direct diagonalization requires knowledge of
all elements of $\rho_{n}$, that is $O\left(d_{s}^{n}\right)$ real numbers.
Projecting $\rho_{n}$ onto its support requires storing and multiplying $2n$
copies of $M$, $n$ copies each of $W$ and $W^{\dagger}$ and 2 copies of
$\Sigma^{-\frac{1}{2}}$ (refer to (54)), that is,
$O\left(2nd_{s}d_{M}^{2}+2d_{M}^{4}\right)=O\left(n\right)$ real numbers.
Therefore, the projection method is exponentially more efficient than direct
diagonalization..

Consider an iuMPS, $\rho_{n}$ given by (7). The support of $\rho_{n}$ is
constructed as follows. Define $G\in\mathbb{C}^{d_{M}^{2}\times d_{M}^{2}}$
by,

$G:=E^{n}=\sum_{s_{1},s_{2},\dots,s_{n}=1}^{d_{s}}\left(M^{s_{n}}\otimes\bar{M}^{s_{n}}\right)\dots{\\\
}\dots\left(M^{s_{2}}\otimes\bar{M}^{s_{2}}\right)\left(M^{s_{1}}\otimes\bar{M}^{s_{1}}\right),$
(53)

which can be numerically computed efficiently. Then, express the matrix
elements of the permutation of $G$, defined by
$\tilde{G}_{ij;i^{\prime}j^{\prime}}:=G_{jj^{\prime};ii^{\prime}}$, in terms
of its spectral decomposition,

$\displaystyle\tilde{G}_{ij;i^{\prime}j^{\prime}}:=\sum_{l,l^{\prime}=1}^{d_{M}^{2}}W_{ij;l}\Sigma_{l;l^{\prime}}\overline{W_{i^{\prime}j^{\prime};l^{\prime}}},$
(54)

which can also be computed efficiently, and where
$i,j,i^{\prime},j^{\prime}\in\left\\{1,2,\dots d_{M}\right\\}$, unitary
$W\in\mathbb{C}^{d_{M}^{2}\times d_{M}^{2}}$, and $\Sigma$ is the diagonal
matrix of eigenvalues of $\tilde{G}$. Equations (53) and (54) are
diagrammatically depicted in Fig. 4. The equation (54) implies that the
support of $\rho_{n}$, denoted as $\mathrm{supp}(\rho_{n})$, is isomorphic to
$\mathcal{H}_{M}\otimes\mathcal{H}_{M}$. We define isometry
$P:\mathcal{H}_{M}\otimes\mathcal{H}_{M}\rightarrow\mathcal{H}^{\otimes
n}_{s}$, with $\mathrm{range}(P)=\mathrm{supp}(\rho_{n})$, defined by its
matrix elements,

$P_{s_{1},s_{2},\dots
s_{n};l}:=\sum_{i,j=1}^{d_{M}}\sum_{k=1}^{d_{M}^{2}}\overline{\left(M^{s_{n}}\dots
M^{s_{2}}M^{s_{1}}\right)_{ij}}\times\\\ \times
W_{ij;k}\Sigma^{-\frac{1}{2}}_{k;l}.$ (55)

The isometry $P$ satisfies

$\displaystyle P^{\dagger}\mathds{1}_{d_{s}^{n}}P=\mathds{1}_{d_{M}^{2}},$
(56)

as diagrammatically verified in Fig. 5. The projection of $\rho_{n}$ onto its
support is given by

$P^{\dagger}\rho_{n}P=\Sigma^{\frac{1}{2}}W^{\dagger}(\sigma\otimes\mathbb{1}_{d_{M}})W\Sigma^{\frac{1}{2}},$
(57)

a $d_{M}^{2}\times d_{M}^{2}$ matrix, with non-zero eigenvalues identical to
those of $\rho_{n}$. Figure 6 depicts $P^{\dagger}\rho_{n}P$. From (57) we
observe that the $d_{M}^{2}\times d_{M}^{2}$ matrix $P^{\dagger}\rho_{n}P$ can
be efficiently computed from $W$ and $\Sigma$, without explicitly computing
$P$, which would be prohibitively expensive. The von Neumann entropy, denoted
by $S\left(\rho_{n}\right)$, is then found by computing the eigenvalues
$\lambda_{i}$ of $P^{\dagger}\rho_{n}P$, and is given by

$\displaystyle
S\left(\rho_{n}\right):=-\sum_{i=1}^{d_{M}^{2}}\lambda_{i}\ln\lambda_{i}.$
(58)

Finally, we repeat the above projection procedure, finding the non-zero
eigenvalues and entropies for each of the four reduced density operators
appearing in the definition of QCMI (10).

Remark:
Alternatively, to calculate the entropy $S(\rho_{R})$ for the reduced density
operator on the continuous region $R$, we could have used the relation (14),
realising that $\rho_{R}$ and

$\displaystyle J^{\dagger}\rho_{R}J$
$\displaystyle=(\mathcal{E}^{|R|}\otimes\mathbb{1}_{M^{\prime}})(|\sqrt{\sigma}\rangle\langle\sqrt{\sigma}|))$
$\displaystyle\in\mathcal{B}(\mathcal{H}_{M}\otimes\mathcal{H}_{M})$

have the same eigenvalues and, thus, entropies. The numerical calculation
would be of comparable complexity to the one we have chosen to employ.

##  5 Numerical verification of Theorem 1

In this section, we present numerical calculations of QCMI for a randomly
generated collection of iuMPS. The main result is given in Fig. 8, which (1)
shows that the QCMI decays approximately exponentially as a function of $|B|$,
and (2) compares its numerical decay rate with the decay rate
$q=2\ln\left(\nu_{\mathrm{gap}}^{-1}\right)$ of the bounding function

$I(A:C|B)\leq Qe^{-q(|B|-K)+2K\ln|B|},$

derived in Theorem 1. Before presenting Fig. 8, we discuss the choice of
parameters for the numerical calculations and the results of benchmarking our
tensor network algorithm against the analytical example in Appendix A.

###  5.1 Choice of parameters

There are 4 parameters in the numerical calculations: sizes of subregions $A$
and $C$ and the dimensions of the physical and virtual spaces $d_{s}$ and
$d_{M}$, respectively. We fix $d_{s}=3$ and $d_{M}=4$, and the sizes of
subregions $A$ and $C$ are held fixed at 1 lattice site each. $\left\lvert
B\right\rvert$ is varied from 2 to 40 in multiples of 2 (to decrease the
complexity of the tensor network algorithm for projecting the reduced density
operators onto their supports).

###  5.2 Benchmarking the numerical calculations

In the example discussed in Appendix A, quantum mutual information (QMI),
defined by

$I(A:C)=S(\rho_{A})+S(\rho_{C})-S(\rho_{AC}),$

is analytically shown to converge to the non-zero finite value,
$I_{\mathrm{th}}:=\lim_{\left\lvert
B\right\rvert\rightarrow\infty}I\left(A:C\right)$, given by,

$\displaystyle I_{\mathrm{th}}$
$\displaystyle=17\left.\ln\left(2\right)\right/16-9\left.\ln\left(3\right)\right/8+5\left.\ln\left(5\right)\right/16,$
$\displaystyle\approx 3.5\times 10^{-3}.$

In Fig. 7, we show the corresponding numerically computed QMI as a function of
$|B|$ in red. We observe that already for $\left\lvert B\right\rvert=26$ the
QMI is identical to $I_{\mathrm{th}}$ down to machine precision
($O(10^{-15})$). This agreement between the analytical and numerical QMI
values benchmarks our tensor network algorithm. With the same algorithm, we
compute the QCMI for the considered example and show the result in black in
Fig. 7. We observe that the QCMI decays exponentially, which is consistent
with Theorem 1.

$\lvert
B\rvert$$\ln\left(I\left(A:C\right)\right)$$\ln\left(I\left(A:C\left\lvert
B\right.\right)\right)$ Figure 7: The numerically computed QMI (in red
circles) and QCMI (in black diamonds) for the example in Appendix A that
converge to a finite value and exponentially to 0, respectively. Left (red)
and right (black) hand side vertical axes represent the natural logarithms of
QMI and QCMI, respectively. The horizontal axes represents the size of the
separating region, $\lvert B\rvert$. The solid light red and gray lines
connecting the red circles and black diamonds guide the readers’ eyes towards
the trends in QMI and QCMI, respectively. The theoretically computed QMI value
for $\left\lvert B\right\rvert\rightarrow\infty$ is denoted by
$I_{\mathrm{th}}$ (dotted black line) and is given in the text. The choice of
parameter values is given in Sec. 5.1.
$|B|-|B_{i}|$$f\left(|B|\right)-f\left(|B_{i}|\right)$ Figure 8: This figure
illustrates the graphs $G_{i}\in\Gamma$ given by (60) as line plots (inset
panel) and as a 2D histogram (main panel). Refer to text in Sec. 5.3 for the
description of the figure.

###  5.3 Main Numerical Result

Let $\mathcal{F}=\\{1,2,...,N\\}$ denote a subset of integers, and take
$i\in\mathcal{F}$ to label each instance of $N$ randomly generated iuMPS
constructed using the methods described in Sec. 4; here $N=6\times 10^{4}$.
Define the function $f_{i}:\mathcal{B}_{i}\rightarrow\mathbb{R},$ where
$i\in\mathcal{F}$ and $\mathcal{B}_{i}:=\\{0,2,4,\cdots,|B_{i}|\\},$ such
that,

$\displaystyle f_{i}(|B|):=\ln I_{i}(A:C|B)/2\ln(\nu_{\mathrm{gap}}^{-1}).$
(59)

Then the inset panel in Fig. 8 displays the collection of graphs
$\Gamma=\\{G_{i}:i\in\mathcal{F}\\}$ where,

$G_{i}:=\\{\left(x-x_{i},f_{i}(x)-f_{i}(x_{i})\right):x=|B|,\\\
x_{i}=|B_{i}|\\}.$ (60)

Note that for graph $G_{i}$ the values of $|B|\leq|B_{i}|$, which is the
maximum size of region $B$. Several graphs could have the same value of
$|B_{i}|$, i.e. the function $i\mapsto|B_{i}|$ is many to one. In practice
$|B_{i}|$ is found by the first instance that the condition
$\exp\left(2\ln\left(\nu_{\mathrm{gap}}^{-1}\right)f_{i}(|B_{i}|)\right)>10^{-k}$
does not hold, where $k=12$ is determined independently by errors in numerical
eigenvalue computations as discussed in Appendix D. In Fig. 8 (in both
panels), we plot blue lines of slope $-1$ to help compare graph data with the
theorem. The origin $\left(0,0\right)$ of both the inset and the main panels
in Fig. 8 corresponds to the maximum dimension $|B|=|B_{i}|$. By comparing the
graphs (black lines) against the blue line of slope $-1$ in the inset panel,
we see that the slopes of most graphs in $\Gamma$ are eventually less than or
equal to $-1$. This indicates, according to Theorem 1, that for most graphs in
$\Gamma$, $|B|=|B_{i}|$ explores the asymptotic regime in the numerical
calculations. In the main panel of Fig. 8, we display the same graph data as a
2D histogram, where the bins of the histogram $H_{i,j}\subset\mathbb{R}^{2}$
are given by,

$H_{i,j}=\left(-\left(2i+1\right),-\left(2i-1\right)\right)\times\left(2j,2j+2\right),$

for $i,j\in\left\\{0,1,2,\dots,19\right\\}$. We observe that for a large
fraction of graphs in $\Gamma$, the slopes are less than or equal to $-1$ and
the occupancy of the bins just below the diagonal of slope $-1$ increases as
$|B|\rightarrow|B_{i}|$.

Normalized decay rate Figure 9:  This figure depicts the cumulative
distribution of the QCMI decay rates normalized by
$2\ln{\nu_{\mathrm{gap}}^{-1}}$ computed (i) for all $G_{i}\in\Gamma$,
illustrated by the dashed black line, and (ii) for subset of $\Gamma$
corresponding to $|B_{i}|=40$, illustrated by the solid black line. The
vertical solid blue line corresponds to the slope of the solid blue line in
Fig. 8. Note that the cumulative distribution extends further than 1.22 on the
horizontal axis, but is not significant to our discussion and is omitted.

In Fig. 9, we further probe the decay of QCMI by plotting the cumulative
distribution of the QCMI decay rates. The decay rates are extracted by linear
regression from the graphs $G_{i}\in\Gamma$. Note that for all the iuMPS in
our generated sample, the constant $K$ is zero, and the logarithmic correction
is absent, which we verify during the computation. This may be observed from
the fact, that in the bound given in Theorem 1, the constant $K$ is non-zero,
only if the Jordan decomposition of the transfer matrix contains non-trivial
Jordan blocks (i.e., contains nilpotent operators) for some eigenvalues with
magnitude less than $1$. (For the eigenvalues of magnitude $1$, Jordan blocks
are trivial, as shown, for example, in the proof of Proposition 3.1 of
[FNW92].) For randomly generated states described in Cases 1-3 of Section 4.1,
non-trivial Jordan blocks will typically be absent, since the eigenvalues (of
magnitude less than 1 for Cases 2 and 3) are typically non-degenerate. Since
we compute the eigenvalues of the transfer matrix, we can check that the
eigenvalues are indeed non-degenerate for all the generated iuMPS.

The dashed black line represents the distribution for all $G_{i}\in\Gamma,$
and the solid black line represents the distribution for graphs $G_{i}$ in the
subset of $\Gamma$ corresponding to $|B_{i}|=\mathrm{max}_{j}|B_{j}|=40$, the
predefined maximum explored in our numerical calculations. Recall that in
other cases $|B_{i}|$ is defined as the value of $|B|$ two units prior to the
value of $|B|$ at which the QCMI reaches the numerical QCMI error. The latter
is determined primarily by errors accrued in calculating the eigenvalues of
the transfer matrix (refer to Appendix D). In other words, the dashed black
line includes those graphs that reach QCMI error for $|B_{i}|<40$, whereas,
the solid black line corresponds to only those graphs that reach QCMI error
for $|B_{i}|\geq 40$. As in Fig. 8, a large fraction of all graphs in $\Gamma$
are consistent with Theorem 1 (dashed black line), because they have decay
rates that obey the asymptotic lower bound. A still higher fraction of graphs
with $|B_{i}|=40$ (solid black line) is consistent with the theorem.
Therefore, with the limited choice of numerical parameters: $\left\lvert
B_{i}\right\rvert\leq 40$ and $(d_{s},d_{M})=(3,4),$ Fig. 8 and Fig. 9 taken
together provide supporting evidence for Theorem 1. In particular Fig. 9
suggests that those example which appear to contradict the theorem reach the
QCMI error for $|B_{i}|<40$, before entering their asymptotic regime, for
which Theorem 1 holds. We suggest that the examples with $|B_{i}|\geq 40$ that
contradict the theorem have also not yet reached their asymptotic regime.

##  6 Conclusion

In this paper we considered a system $ABC$ as shown in Fig.2(a) and (b),
prepared in an infinite uniform matrix product state with density operator
$\rho_{ABC}$, illustrated in Fig.2(c). We proved that the QCMI, denoted
$I(A:C|B)$, is bounded by an exponentially decaying function in the asymptotic
limit of large size $|B|$ of the separating region. By taking advantage of the
more direct methods of estimation than in [SK22], we have improved the
exponent of the asymptotic bounding function for QCMI by a factor of 4. We
have also shown that $\rho_{ABC}$ is close to an exact quantum Markov chain
$\tilde{\rho}_{ABC}$ in the trace distance, which is also bounded by an
exponentially decaying function. Compared to [SK22], the exponent of the
latter function is improved by a factor of 2. We speculate that the achieved
asymptotic bounds are tight. The theory was supported by numerical experiments
in which the QCMI is computed for a collection of matrix product states
generated by sampling their generating isometries with respect to Haar
measure.

## Acknowledgement

We are grateful to Glen Evenbly for useful discussions on the numerical
methods involving tensor networks.

## References

  * [AF04] R Alicki and M Fannes  “Continuity of quantum conditional information”  IOP Publishing, 2004, pp. L55–L57  DOI: 10.1088/0305-4470/37/5/l01
  * [BH15] Fernando G… Brandão and Michał Horodecki  “Exponential Decay of Correlations Implies Area Law”  In _Communications in Mathematical Physics_ 333.2, 2015, pp. 761–798  DOI: 10.1007/s00220-014-2213-8
  * [BK19] Fernando G… Brandão and Michael J. Kastoryano  “Finite Correlation Length Implies Efficient Preparation of Quantum Thermal States”  In _Communications in Mathematical Physics_ 365.1, 2019, pp. 1–16  DOI: 10.1007/s00220-018-3150-8
  * [Bru+10] Wojciech Bruzda et al.  “Universality of spectra for interacting quantum chaotic systems”  In _Phys. Rev. E_ 81 American Physical Society, 2010, pp. 066209  DOI: 10.1103/PhysRevE.81.066209
  * [Cho75] Man-Duen Choi  “Completely positive linear maps on complex matrices”  In _Linear Algebra and its Applications_ 10.3, 1975, pp. 285–290  DOI: https://doi.org/10.1016/0024-3795(75)90075-0
  * [CKB21] Chi-Fang Chen, Kohtaro Kato and Fernando G… Brandão  “Matrix Product Density Operators: when do they have a local parent Hamiltonian?”, 2021  arXiv:2010.14682 [quant-ph]
  * [Fla+17] Steven T. Flammia, Jeongwan Haah, Michael J. Kastoryano and Isaac H. Kim  “Limits on the storage of quantum information in a volume of space”  In _Quantum_ 1 Verein zur Förderung des Open Access Publizierens in den Quantenwissenschaften, 2017, pp. 4  DOI: 10.22331/q-2017-04-25-4
  * [FNW92] M. Fannes, B. Nachtergaele and R.. Werner  “Finitely correlated states on quantum spin chains”  In _Communications in Mathematical Physics_ 144.3, 1992, pp. 443–490  DOI: 10.1007/BF02099178
  * [FNW94] M. Fannes, B. Nachtergaele and R.. Werner  “Finitely Correlated Pure States”  In _Journal of Functional Analysis_ 120.2, 1994, pp. 511–534  DOI: https://doi.org/10.1006/jfan.1994.1041
  * [FR15] Omar Fawzi and Renato Renner  “Quantum Conditional Mutual Information and Approximate Markov Chains”  In _Communications in Mathematical Physics_ 340.2, 2015, pp. 575–611  DOI: 10.1007/s00220-015-2466-x
  * [GJN18] Carlos E. González-Guillén, Marius Junge and Ion Nechita  “On the spectral gap of random quantum channels”  arXiv, 2018  DOI: 10.48550/ARXIV.1811.08847
  * [Gut+20] Jiří Guth Jarkovský, András Molnár, Norbert Schuch and J. Cirac  “Efficient Description of Many-Body Systems with Matrix Product Density Operators”  In _PRX Quantum_ 1.1 American Physical Society, 2020, pp. 010304  DOI: 10.1103/PRXQuantum.1.010304
  * [Has07] M.. Hastings  “Random unitaries give quantum expanders”  In _Phys. Rev. A_ 76 American Physical Society, 2007, pp. 032315  DOI: 10.1103/PhysRevA.76.032315
  * [Hay+04] Patrick Hayden, Richard Jozsa, Dénes Petz and Andreas Winter  “Structure of States Which Satisfy Strong Subadditivity of Quantum Entropy with Equality”  In _Communications in Mathematical Physics_ 246.2, 2004, pp. 359–374  DOI: 10.1007/s00220-004-1049-z
  * [HP20] Patrick Hayden and Geoffrey Penington  “Approximate Quantum Error Correction Revisited: Introducing the Alpha-Bit”  In _Communications in Mathematical Physics_ 374.2, 2020, pp. 369–432  DOI: 10.1007/s00220-020-03689-1
  * [Jun+18] Marius Junge et al.  “Universal Recovery Maps and Approximate Sufficiency of Quantum Relative Entropy”  In _Annales Henri Poincaré_ 19.10, 2018, pp. 2955–2978  DOI: 10.1007/s00023-018-0716-0
  * [KAA21] Tomotaka Kuwahara, Álvaro M. Alhambra and Anurag Anshu  “Improved Thermal Area Law and Quasilinear Time Algorithm for Quantum Gibbs States”  In _Phys. Rev. X_ 11 American Physical Society, 2021, pp. 011047  DOI: 10.1103/PhysRevX.11.011047
  * [KB19] Kohtaro Kato and Fernando G… Brandão  “Quantum Approximate Markov Chains are Thermal”  In _Communications in Mathematical Physics_ 370.1, 2019, pp. 117–149  DOI: 10.1007/s00220-019-03485-6
  * [Kim17] Isaac Kim  “Markovian Matrix Product Density Operators : Efficient computation of global entropy”, 2017  arXiv:1709.07828 [quant-ph]
  * [KKB20] Tomotaka Kuwahara, Kohtaro Kato and Fernando G… Brandão  “Clustering of Conditional Mutual Information for Quantum Gibbs States above a Threshold Temperature”  In _Phys. Rev. Lett._ 124 American Physical Society, 2020, pp. 220601  DOI: 10.1103/PhysRevLett.124.220601
  * [KSW08] Dennis Kretschmann, Dirk Schlingemann and Reinhard F. Werner  “The Information-Disturbance Tradeoff and the Continuity of Stinespring’s Representation”  In _IEEE Transactions on Information Theory_ 54.4, 2008, pp. 1708–1717  DOI: 10.1109/TIT.2008.917696
  * [KSZ93] A Klümper, A Schadschneider and J Zittartz  “Matrix Product Ground States for One-Dimensional Spin-1 Quantum Antiferromagnets”  IOP Publishing, 1993, pp. 293–297  DOI: 10.1209/0295-5075/24/4/010
  * [Kuk+21] Ryszard Kukulski et al.  “Generating random quantum channels”  In _Journal of Mathematical Physics_ 62.6, 2021, pp. 062201  DOI: 10.1063/5.0038838
  * [LP21] Cécilia Lancien and David Pérez-García  “Correlation Length in Random MPS and PEPS”  In _Annales Henri Poincaré_ 23.1 Springer ScienceBusiness Media LLC, 2021, pp. 141–222  DOI: 10.1007/s00023-021-01087-4
  * [LR68] Oscar E. Lanford and Derek W. Robinson  “Mean Entropy of States in Quantum-Statistical Mechanics”  In _Journal of Mathematical Physics_ 9.7, 1968, pp. 1120–1125  DOI: 10.1063/1.1664685
  * [LR73] Elliott H. Lieb and Mary Beth Ruskai  “A Fundamental Property of Quantum-Mechanical Entropy”  In _Phys. Rev. Lett._ 30 American Physical Society, 1973, pp. 434–436  DOI: 10.1103/PhysRevLett.30.434
  * [LR73a] Elliott H. Lieb and Mary Beth Ruskai  “Proof of the strong subadditivity of quantum-mechanical entropy”  In _Journal of Mathematical Physics_ 14.12, 1973, pp. 1938–1941  DOI: 10.1063/1.1666274
  * [Pér+07] D. Pérez-García, F. Verstraete, M.. Wolf and J.. Cirac  “Matrix Product State Representations”  In _Quantum Info. Comput._ 7.5 Paramus, NJ: Rinton Press, Incorporated, 2007, pp. 401–430
  * [Pet86] Dénes Petz  “Sufficient subalgebras and the relative entropy of states of a von Neumann algebra”  In _Communications in Mathematical Physics_ 105.1, 1986, pp. 123–131  DOI: 10.1007/BF01212345
  * [Pet88] Dénes Petz  “Sufficiency of channels over von Neumann algebras”  In _The Quarterly Journal of Mathematics_ 39.1, 1988, pp. 97–108  DOI: 10.1093/qmath/39.1.97
  * [Rus02] Mary Beth Ruskai  “Inequalities for quantum entropy: A review with conditions for equality”  In _Journal of Mathematical Physics_ 43.9, 2002, pp. 4358–4375  DOI: 10.1063/1.1497701
  * [Sch+07] C. Schön et al.  “Sequential generation of matrix-product states in cavity QED”  In _Phys. Rev. A_ 75 American Physical Society, 2007, pp. 032311  DOI: 10.1103/PhysRevA.75.032311
  * [SFR16] David Sutter, Omar Fawzi and Renato Renner  “Universal recovery map for approximate Markov chains”  In _Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences_ 472.2186, 2016, pp. 20150623  DOI: 10.1098/rspa.2015.0623
  * [SK22] Pavel Svetlichnyy and T… Kennedy  “Decay of quantum conditional mutual information for purely generated finitely correlated states”  In _Journal of Mathematical Physics_ 63.7, 2022, pp. 072201  DOI: 10.1063/5.0085358
  * [SM16] Brian Swingle and John McGreevy  “Mixed $s$-sourcery: Building many-body states using bubbles of nothing”  In _Phys. Rev. B_ 94 American Physical Society, 2016, pp. 155125  DOI: 10.1103/PhysRevB.94.155125
  * [SRW15] Oleg Szehr, David Reeb and Michael M. Wolf  “Spectral Convergence Bounds for Classical and Quantum Markov Processes”  In _Communications in Mathematical Physics_ 333.2, 2015, pp. 565–595  DOI: 10.1007/s00220-014-2188-5
  * [Sti55] W. Stinespring  “Positive Functions on C*-Algebras”  In _Proceedings of the American Mathematical Society_ 6.2 American Mathematical Society, 1955, pp. 211–216  URL: http://www.jstor.org/stable/2032342
  * [Sut18] David Sutter  “Approximate Quantum Markov Chains”  Springer Cham, 2018  DOI: 10.1007/978-3-319-78732-9
  * [VGC04] F. Verstraete, J.. García-Ripoll and J.. Cirac  “Matrix Product Density Operators: Simulation of Finite-Temperature and Dissipative Systems”  In _Phys. Rev. Lett._ 93 American Physical Society, 2004, pp. 207204  DOI: 10.1103/PhysRevLett.93.207204
  * [Vid03] Guifré Vidal  “Efficient Classical Simulation of Slightly Entangled Quantum Computations”  In _Phys. Rev. Lett._ 91 American Physical Society, 2003, pp. 147902  DOI: 10.1103/PhysRevLett.91.147902
  * [VMC08] F. Verstraete, V. Murg and J.I. Cirac  “Matrix product states, projected entangled pair states, and variational renormalization group methods for quantum spin systems”  In _Advances in Physics_ 57.2 Taylor & Francis, 2008, pp. 143–224  DOI: 10.1080/14789940801912366
  * [Whi92] Steven R. White  “Density matrix formulation for quantum renormalization groups”  In _Phys. Rev. Lett._ 69 American Physical Society, 1992, pp. 2863–2866  DOI: 10.1103/PhysRevLett.69.2863
  * [Wil17] Mark M. Wilde  “Quantum Information Theory”  Cambridge University Press, 2017  DOI: 10.1017/9781316809976
  * [Win16] Andreas Winter  “Tight Uniform Continuity Bounds for Quantum Entropies: Conditional Entropy, Relative Entropy Distance and Energy Constraints”  In _Communications in Mathematical Physics_ 347.1, 2016, pp. 291–313  DOI: 10.1007/s00220-016-2609-8
  * [Wol+08] Michael M. Wolf, Frank Verstraete, Matthew B. Hastings and J. Cirac  “Area Laws in Quantum Systems: Mutual Information and Correlations”  In _Phys. Rev. Lett._ 100 American Physical Society, 2008, pp. 070502  DOI: 10.1103/PhysRevLett.100.070502
  * [ZCL21] Leo Zhou, Soonwon Choi and Mikhail D. Lukin  “Symmetry-protected dissipative preparation of matrix product states”  In _Phys. Rev. A_ 104 American Physical Society, 2021, pp. 032418  DOI: 10.1103/PhysRevA.104.032418

##  Appendix A An example: convergence of QMI to a non-zero value

In this appendix we specify the parameters of the iuMPS used in Fig.7 and
Fig.8, and the asymptotic limit of its QMI

$I(A:C)=S(\rho_{A})+S(\rho_{C})-S(\rho_{AC}),$

which we use to benchmark our numerical computations.

The iuMPS is generated by

$M^{s}=\begin{bmatrix}M_{1}^{s}&0\\\ 0&M_{2}^{s}\end{bmatrix},$

i.e., belongs to Case 2 defined in Sec. 4.1, with submatrices specified by,

$\displaystyle M_{1}^{1}=M_{2}^{3}=\begin{bmatrix}0&-\frac{1}{\sqrt{2}}\\\
0&0\end{bmatrix},$ $\displaystyle
M_{1}^{2}=M_{2}^{2}=\begin{bmatrix}-\frac{1}{\sqrt{2}}&0\\\
0&\frac{1}{\sqrt{2}}\end{bmatrix},$ (61) $\displaystyle
M_{1}^{3}=M_{2}^{1}=\begin{bmatrix}0&0\\\ \frac{1}{\sqrt{2}}&0\end{bmatrix},$

and by the fixed point $|\sigma\rangle=(1/2,1/2,1/2,1/2)^{T}$.

As was shown in Appendix D of [SK22], if we set $|A|=|C|=1$, then, when
$|B|\rightarrow\infty$, the reduced density operators $\rho_{AC}$, $\rho_{A}$,
and $\rho_{C}$ converge to

$\displaystyle\rho_{AC}$ $\displaystyle=\frac{1}{32}\begin{pmatrix}1&0&0\\\
0&2&0\\\ 0&0&1\end{pmatrix}\otimes\begin{pmatrix}1&0&0\\\ 0&2&0\\\
0&0&1\end{pmatrix}$ (62) $\displaystyle+\frac{1}{32}\begin{pmatrix}1&0&0\\\
0&1&0\\\ 0&0&2\end{pmatrix}\otimes\begin{pmatrix}1&0&0\\\ 0&1&0\\\
0&0&2\end{pmatrix},$ $\displaystyle\rho_{A}=$
$\displaystyle\frac{1}{8}\begin{pmatrix}2&0&0\\\ 0&3&0\\\
0&0&3\end{pmatrix},\;\rho_{C}=\frac{1}{8}\begin{pmatrix}2&0&0\\\ 0&3&0\\\
0&0&3\end{pmatrix}.$

By direct calculation, the QMI $I(A:C)$, converges to the value
$I_{\mathrm{th}}=17\left.\ln\left(2\right)\right/16-9\left.\ln\left(3\right)\right/8+5\left.\ln\left(5\right)\right/16\approx
3.5\times 10^{-3}$. In contrast, as is proven by Theorem 1 or by Theorem II.1
of [SK22], the QCMI converges to zero.

##  Appendix B Taylor expansion with remainder

In this appendix we provide details of the derivation of the bound (19), that
we have omitted in Subsection 3.1.

Recall that we are aiming to bound the QCMI

$\displaystyle
I(A:C|B)=\mathrm{tr}_{M}\left(\rho_{ABC}\left(\ln\rho_{ABC}+\ln\rho_{B}-\ln\rho_{AB}-\ln\rho_{BC}\right)\right)$
(63)

in terms of $\|\Delta\rho_{ABC}\|_{1}$, where
$\Delta\rho_{ABC}=\rho_{ABC}-\tilde{\rho}_{ABC}$. We achieve this by making
the Taylor expansion of all the logarithms, which, with an explicit remainder
term, takes the form,

$\displaystyle\ln(\tilde{\rho}+\Delta\rho)=\ln\tilde{\rho}$
$\displaystyle+\int_{0}^{+\infty}ds\left(\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}-\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+s\mathbb{1}}\right.$
(64)
$\displaystyle\quad+\left.\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\right),$

where $-1<t<1$, and where we dropped index that can be any of $ABC$, $AB$,
$BC$, or $B$.

We have already shown in Subsection 3.1, that the contributions to $I(A:C|B)$
of the zeroth and the first order in $\|\Delta\rho_{ABC}\|_{1}$ vanish. We
have also estimated the contribution of any of the four terms of the second
order in (13).

Now we estimate the contribution of the remainder term in the expansion (64),

$\displaystyle\left|\mathrm{tr}\left((\tilde{\rho}+\Delta\rho)\int_{0}^{+\infty}ds\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\Delta\rho\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\right)\right|$
(65)
$\displaystyle\qquad\qquad\qquad\leq\int_{0}^{+\infty}ds\left\|\frac{1}{\tilde{\rho}+t\Delta\rho+s\mathbb{1}}\right\|^{4}\|\Delta\rho\|_{1}^{3}$
$\displaystyle\qquad\qquad\qquad\leq\int_{0}^{+\infty}\frac{ds}{(\lambda^{\prime}_{\mathrm{min}}+s)^{4}}\|\Delta\rho\|_{1}^{3}$
$\displaystyle\qquad\qquad\qquad=\frac{\|\Delta\rho\|^{3}_{1}}{3\lambda^{\prime
3}_{\mathrm{min}}},$

where $\lambda^{\prime}_{\mathrm{min}}$ is the smallest eigenvalue of
$\tilde{\rho}+t\Delta\rho$.

From (18) in Subsection 3.1 and (44) in Subsection 3.2 it follows that for an
injective iuMPS,

$\displaystyle|\lambda_{\mathrm{min}}-\sigma^{2}_{\mathrm{min}}|\leq\|\rho_{ABC}-\rho^{\prime}_{ABC}\|+\|\Delta\rho\|\leq
c_{2}d_{M}\left(1+\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}\right)|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}.$
(66)

A similar calculation for $\tilde{\rho}_{R}+t\Delta\rho_{R}$ leads to the
bound

$\displaystyle|\lambda^{\prime}_{\mathrm{min}}-\sigma^{2}_{\mathrm{min}}|\leq\|\rho_{ABC}-\rho^{\prime}_{ABC}\|+2\|\Delta\rho\|\leq
c_{2}d_{M}\left(1+2\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}\right)|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}.$
(67)

Combining (13), (65), (66) and (67), and recalling that there are four
logarithms contributing to $I(A:C|B)$ (63), we estimate,

$\displaystyle I(A:C|B)\leq
4\left(\frac{\|\Delta\rho_{ABC}\|^{2}_{1}}{2\lambda_{\mathrm{min}}}+\frac{\|\Delta\rho_{ABC}\|^{3}_{1}}{3\lambda^{\prime
3}_{\mathrm{min}}}\right)=\frac{2\|\Delta\rho_{ABC}\|^{2}_{1}}{\lambda_{\mathrm{min}}}\left(1+\frac{2\lambda_{\mathrm{min}}}{3\lambda^{\prime
3}_{\mathrm{min}}}\|\Delta\rho_{ABC}\|_{1}\right),$

where we have redefined $\lambda_{\mathrm{min}}$ to be the smallest eigenvalue
of all $\tilde{\rho}_{R}$, where $R$ can be any of $ABC$, $AB$, $BC$, or $B$.
We have similarly redefined $\lambda^{\prime}_{\mathrm{min}}$, with
$\tilde{\rho}_{R}+t\Delta\rho_{R}$ replacing $\tilde{\rho}_{R}$. We have also
used the inequality $\|\Delta\rho_{R}\|_{1}\leq\|\Delta\rho_{ABC}\|_{1}$.

From (66) and (67) we observe that $\lambda_{\mathrm{min}}$ and
$\lambda^{\prime}_{\mathrm{min}}$ both converge to
$\sigma^{2}_{\mathrm{min}}$, while (from (44)) the trace distance
$\|\Delta\rho_{ABC}\|_{1}$ vanishes as $|B|\rightarrow\infty$. Then, for large
enough $|B|$,

$\frac{2\lambda_{\mathrm{min}}}{3\lambda^{\prime
3}_{\mathrm{min}}}\|\Delta\rho_{ABC}\|_{1}\leq 1.$ (68)

This leads to the bound

$I(A:C|B)\leq\frac{4\|\Delta\rho_{ABC}\|^{2}_{1}}{\lambda_{\mathrm{min}}}.$

For large enough $|B|$, for an injective iuMPS, from (66) we observe that we
can guarantee the bound

$\lambda_{\mathrm{min}}\geq\frac{\sigma^{2}_{\mathrm{min}}}{2}.$ (69)

Then,

$I(A:C|B)\leq\frac{8}{\sigma^{2}_{\mathrm{min}}}\|\Delta\rho_{ABC}\|^{2}_{1}.$

##  Appendix C Sufficient conditions on $|B|$

In this appendix we combine together all the requirements imposed on the size
$|B|$ of the separating region in order to get a unified sufficient condition.

First, we require $|B|$ to be large enough, that the bound (9) is satisfied.
We derive such a condition in Appendix F,

$\begin{cases}|B|>2K,\\\
\left(\frac{|B|}{K}\right)^{K}\nu_{\mathrm{gap}}^{|B|-K}\geq(K_{\nu_{i}}+1)\left(\frac{|B|e}{K_{\nu_{i}}}\right)^{K_{\nu_{i}}}|\nu_{i}|^{|B|-K_{\nu_{i}}},\;\mathrm{for\;all\;}\nu_{i}\in\mathcal{O}\backslash\mathcal{O}_{*}.\end{cases}$

The constants $K$ and $K_{\nu_{i}}$ may be extracted from the Jordan
decomposition of $E$, as described in Appendix F.

Second, we consider the requirement imposed on $|B|$ in Appendix B in (69).
According to (66), this requirement may be satisfied by setting

$c_{2}d_{M}\left(1+\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}\right)|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{2}\sigma^{2}_{\mathrm{min}}.$

It will be more convenient for us to satisfy a stronger condition

$c_{2}d_{M}\left(1+2\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}\right)|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{2}\sigma^{2}_{\mathrm{min}},$
(70)

which may be guaranteed by satisfying an even stronger condition

$3c_{2}d_{M}\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{2}\sigma^{2}_{\mathrm{min}},$

where we have used $\sqrt{2d_{M}/\sigma_{\mathrm{min}}}>1$. This results in
the condition

$|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{6\sqrt{2}}\frac{\sigma^{\frac{5}{2}}_{\mathrm{min}}}{c_{2}d_{M}^{\frac{3}{2}}}.$
(71)

From the condition (70), it follows that

$|\lambda^{\prime}_{\mathrm{min}}-\sigma^{2}_{\mathrm{min}}|\leq\frac{1}{2}\sigma^{2}_{\mathrm{min}},$

which we employ in the sequel.

Now we consider the third requirement imposed on $|B|$ in Appendix B, the
equation (68), which leads to the condition

$c_{2}d_{M}\sqrt{\frac{2d_{M}}{\sigma_{\mathrm{min}}}}|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq
3\left(\frac{3}{2}\sigma^{2}_{\mathrm{min}}\right)^{3}\sigma^{-2}_{\mathrm{min}}.$

Expressing the above inequality in the same form as (71), we arrive at the
formula

$|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{81}{8\sqrt{2}}\frac{\sigma^{\frac{9}{2}}_{\mathrm{min}}}{c_{2}d^{\frac{3}{2}}_{M}}.$
(72)

We imposed a further condition on $|B|$, namely $\|\Delta C^{\prime}\|\leq 1$,
in order to derive (41), with $\Delta C^{\prime}$ defined by (30) and (38). It
was superseded by the stronger condition $\|\Delta C^{\prime}\|\leq 2^{-3/2}$,
which we used to guarantee the inequality (42). By using (31), this condition
leads to the inequality

$\frac{1}{\sigma_{\mathrm{min}}}c_{2}d_{M}|B|^{K}\nu^{|B|}_{\mathrm{gap}}\leq\frac{1}{2\sqrt{2}}.$

Expressing this inequality in the same form as (71), we arrive at the
condition

$|B|^{K}\nu^{|B|}_{\mathrm{gap}}\leq\frac{1}{2\sqrt{2}}\frac{\sigma_{\mathrm{min}}}{c_{2}d_{M}}.$

This condition is, however, redundant, since it is strictly weaker than (71).

Thus, combining (71) and (72) together, we obtain

$|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{6\sqrt{2}}\frac{\sigma^{\frac{5}{2}}_{\mathrm{min}}}{c_{2}d_{M}^{\frac{3}{2}}}\min\left\\{1,\frac{243}{4}\sigma^{2}_{\mathrm{min}}\right\\}.$

We have imposed one more condition on $|B|$, in Subsection 3.1, by requiring
$\dim\mathcal{H}_{B}\geq d^{2}_{M}$,

$|B|\geq 2\frac{\ln d_{M}}{\ln d_{s}}.$

Thus, the size $|B|$ should be large enough to satisfy the conditions

$\begin{cases}|B|>2K,\\\
\left(\frac{|B|}{K}\right)^{K}\nu_{\mathrm{gap}}^{|B|-K}\geq(K_{\nu_{i}}+1)\left(\frac{|B|e}{K_{\nu_{i}}}\right)^{K_{\nu_{i}}}|\nu_{i}|^{|B|-K_{\nu_{i}}},\;\mathrm{for\;all\;}\nu_{i}\in\mathcal{O}\backslash\mathcal{O}_{*},\\\
|B|^{K}\nu_{\mathrm{gap}}^{|B|-K}\leq\frac{1}{6\sqrt{2}}\frac{\sigma^{\frac{5}{2}}_{\mathrm{min}}}{c_{2}d_{M}^{\frac{3}{2}}}\min\left\\{1,\frac{243}{4}\sigma^{2}_{\mathrm{min}}\right\\},\\\
|B|\geq 2\frac{\ln d_{M}}{\ln d_{s}}.\end{cases}$ (73)

The constant $c_{2}$ is specified in Appendix F.

We remind the reader that the constraints (73) should not be taken as
necessary to prove Theorem 1, and are only sufficient (and convenient) to
guarantee validity of the steps in the proof.

##  Appendix D An example of errors in eigenvalue calculations

Errors resulting from the numerical computation of eigenvalues of reduced
density operators dictate the smallest value of QCMI that we can reliably
compute. In turn, this restricts the maximum size of the separating region
$|B|$ we can explore. Since
$I(A:C|B)=S(\rho_{AB})+S(\rho_{BC})-S(\rho_{ABC})-S(\rho_{B})$, the error
incurred in calculating QCMI is estimated as

$\delta I(A:C|B)\sim\delta S,$

where,

$\delta S:=\mathrm{max}\left\\{\delta S\left(\rho_{ABC}\right),\delta
S\left(\rho_{B}\right),\delta S\left(\rho_{AB}\right),\delta
S\left(\rho_{BC}\right)\right\\},$

is the maximum error in a calculation of the entropies. We cannot reliably
compute values of QCMI smaller than $\delta I$, and thus are forced to stop
calculations when the size of the separating region $|B|$ becomes large enough
that $I(A:C|B)$ decays to values of the order $O(\delta I)$.

As discussed in Section 4.2, each of the density operators $\rho_{ABC}$,
$\rho_{AB}$, $\rho_{BC}$, and $\rho_{B}$ is supported on a subspace of
dimension $d_{M}^{2}$ so

$\displaystyle S=-\sum_{i=1}^{d_{M}^{2}}\lambda_{i}\ln\lambda_{i}.$

Hence, we may estimate the error in calculating entropies as

$\displaystyle\delta S$
$\displaystyle\leq\sum_{i=1}^{d_{M}^{2}}\delta\lambda\ln\lambda_{i}^{-1}-\delta\lambda\lesssim\sum_{i=1}^{d_{M}^{2}}\delta\lambda\ln\lambda_{i}^{-1}\lesssim
d_{M}^{2}\delta\lambda\ln\lambda_{\mathrm{min}}^{-1},$

where $\lambda_{\mathrm{min}}:=\min_{i}\lambda_{i}$, and $\delta\lambda>0$ is
the maximum error in computing the eigenvalues. In practice, we use the
MATLAB’s built in eigenvalue calculation function “eig” to compute
eigenvalues, with an error of order $O\left(10^{-14}\right),$ as explained in
the following paragraph. Hence we find

$\delta I\sim\delta S\sim 4^{2}\cdot 10^{-14}\cdot\ln 10^{14}\sim 10^{-12}.$
(74)

In our numerical calculations, refering to (59) and below (60), for each
$i\in\mathcal{F}$, we denote $\left\lvert B_{i}\right\rvert$ to be the maximum
size of the separating region such that
$\exp\left(2\ln\left(\nu_{\mathrm{gap}}^{-1}\right)f_{i}(|B_{i}|)\right)>10^{-k}.$
From the above discussion, we can reliably compute QCMI as small as
$O\left(10^{-12}\right)$, and, thus, we take $k=12$ for our numerical
calculations.

To gauge the size of errors incurred in eigenvalue calculations we compare
numerically computed values against an exact result. We consider the transfer
matrix eigenvalue problem of Case 1 in Sec. 4.1 and recall that one of the
eigenvalues $\nu_{1}=1.$ The numerically computed value incurs an error of
$O(10^{-14})$. This estimate of the error is derived empirically from Fig.
10(a), which depicts the absolute difference between $\nu_{1}$, and its
numerically computed value for $2\times 10^{4}$ different instances of
transfer matrices.

##  Appendix E Spectrum of the transfer matrix

In this appendix we discuss the spectrum for the transfer matrix of an iuMPS,
which is sampled with respect to the Haar measure (as defined in Section 4).

For the purpose of comparing theoretically expected eigenvalues of the
transfer matrix to numerically computed eigenvalues, both denoted by $\nu_{i}$
and delineated by context in what follows, we assume the ordering
$|\nu_{1}|\geq|\nu_{2}|\geq\dots\geq|\nu_{d_{M}^{2}}|$. Recall that we define

$\nu_{\mathrm{gap}}=\max_{|\nu_{i}|<1}|\nu_{i}|.$

In Fig. 10(a) we depict the magnitude of difference between 1, the known
largest eigenvalue, and $|\nu_{1}|$, the largest numerically computed value.
As mentioned in Appendix D, we infer from the figure that eigenvalue
computations incur a numerical error of the order $O(10^{-14})$.

In Fig. 10(b), the magnitude of the difference between $|\nu_{1}|$ and
$|\nu_{2}|$ is shown. We observe a gap between the computed eigenvalue
magnitudes that clearly exceeds the numerical error. The gap is, however, an
artifact of the finite sampling of Haar random isometries used to generate the
iuMPS’s, as a simple example illustrates. Consider the family of quantum
channels $\mathcal{E}:\mathbb{C}^{4\times 4}\rightarrow\mathbb{C}^{4\times
4}$, generated by the dilating isometry
$V:\mathbb{C}^{4}\rightarrow\mathbb{C}^{3}\otimes\mathbb{C}^{4}$, defined by

$\displaystyle V=$
$\displaystyle\sqrt{\beta}|1\rangle\otimes|-\rangle\langle+|\otimes\mathbb{1}_{2}+\sqrt{1-\beta}|0\rangle\otimes\mathbb{1}_{2}\otimes\mathbb{1}_{2}-\sqrt{\beta}|-1\rangle\otimes|+\rangle\langle-|\otimes\mathbb{1}_{2},$

where we have decomposed
$\mathbb{C}^{4}\cong\mathbb{C}^{2}\otimes\mathbb{C}^{2}$, the bases of
$\mathbb{C}^{2}$ and $\mathbb{C}^{3}$ are denoted by
$\\{|+\rangle,|-\rangle\\}$ and $\\{|1\rangle,|0\rangle,|-1\rangle\\}$,
respectively, and where $0\leq\beta\leq 1$. The isometry induces the quantum
channel

$\mathcal{E}=\mathrm{tr}_{3}(VXV^{\dagger})=(1-\beta)\,\mathrm{id}_{2}\otimes\mathrm{id}_{2}+\beta\,\mathcal{E}_{1}\otimes\mathrm{id}_{2},$

where $\mathrm{id}_{2}:\mathbb{C}^{2\times 2}\rightarrow\mathbb{C}^{2\times
2}$ is the identity channel, and the quantum channel
$\mathcal{E}_{1}:\mathbb{C}^{2\times 2}\rightarrow\mathbb{C}^{2\times 2}$ is
defined by

$\displaystyle\mathcal{E}_{1}(|+\rangle\langle+|)$
$\displaystyle=|-\rangle\langle-|,$
$\displaystyle\mathcal{E}_{1}(|-\rangle\langle-|)$
$\displaystyle=-|+\rangle\langle+|,$
$\displaystyle\mathcal{E}_{1}(|+\rangle\langle-|)$
$\displaystyle=\mathcal{E}_{1}(|-\rangle\langle+|)=0.$

For this quantum channel, we may readily compute

$|\nu_{1}|-|\nu_{2}|=2\beta,$

which may be arbitrarily close to $0$, thus closing the gap in Fig. 10(b).
These channels belong to the family represented by _Case 1_ in Section 4.1.
Our numerical results, however, suggest that the probability of sampling
channels of this kind is small. The spectrum of random quantum channels (Haar-
distributed or otherwise) has been studied in previous works [GJN18, LP21,
Kuk+21, Bru+10, Has07], including the gap in the spectrum of quantum channels.
The presence of a gap in a typical - with a particular sense of this word
depending on the theoretical study - quantum channel was proven, for example,
for Gaussian distributed quantum channels [LP21], and Haar-distributed quantum
channels [GJN18], in an asymptotic limit of large dimensions. These settings
are closely related to this paper, but our numerical computations involve
iuMPS with the virtual dimension $d_{M}=4$ and physical dimension $d_{s}=3$,
both of which are too small for the results of [GJN18] to be applicable.
Nevertheless, our numerical experiments indicate that the gap in the spectrum
of Haar-random quantum channels is present even for such small dimensions.

In Fig. 10(c), we show the magnitude of the numerically computed difference
between $|\nu_{2}|$ and $|\nu_{3}|$. We observe that, in contrast to the
difference $||\nu_{1}|-|\nu_{2}||$ in Fig. 10(b), for a non-zero fraction of
the sampled transfer matrices the second and the third largest eigenvalue
magnitudes are equal (or differ by less than numerical error). We speculate
that these occurrences might be explained as follows. Recall that eigenvalues
of a transfer matrix come in complex conjugate pairs. Hence, if we sample a
transfer matrix with complex eigenvalues, we will necessarily find complex
conjugate pairs of equal modulus. Thus, together, Figs. 10(b) and (c) indicate
that, using the Haar measure, it is possible to sample with non-vanishing
probability transfer matrices that have complex eigenvalues of magnitude
strictly less than 1 (Fig. 10(c)); the probability of sampling transfer
matrices with more than one eigenvalue of magnitude one (including cases of
non-real eigenvalues) is, however, vanishingly small (Fig. 10(b)). (See
[Bru+10] for a relevant study of distributions of eigenvalues of random
quantum channels.)

The gap in the sampled values of $||\nu_{2}|-|\nu_{3}||$ is, as in the case of
$||\nu_{1}|-|\nu_{2}||$, a numerical artifact. For the collection of transfer
matrices generated by the one-parameter family of isometries

$\displaystyle V=$
$\displaystyle\left((1-\beta)\sqrt{\frac{2}{3}}+\beta\frac{\sqrt{3}}{2}\right)|1\rangle\otimes|-\rangle\langle+|\otimes\mathbb{1}_{2}$
$\displaystyle-\left((1-\beta)\frac{1}{\sqrt{3}}+\beta\frac{1}{2}\right)|0\rangle\otimes|+\rangle\langle+|\otimes\mathbb{1}_{2}$
$\displaystyle+\left((1-\beta)\frac{1}{\sqrt{3}}+\beta\frac{1}{2}\right)|0\rangle\otimes|-\rangle\langle-|\otimes\mathbb{1}_{2}$
$\displaystyle-\left((1-\beta)\sqrt{\frac{2}{3}}+\beta\frac{\sqrt{3}}{2}\right)|-1\rangle\otimes|+\rangle\langle-|\otimes\mathbb{1}_{2},$

with $0\leq\beta\leq 1$, the analytically calculated value,

$|\nu_{2}|-|\nu_{3}|=\frac{\sqrt{6}-2}{\sqrt{3}}\beta+O(\beta^{2}),$

may be arbitrarily close to zero.

Instance$\left\lvert
1-\left\lvert\nu_{1}\right\rvert\right\rvert$$\left(\mathrm{a}\right)$

Instance$\left\lvert\left\lvert\nu_{1}\right\rvert-\left\lvert\nu_{2}\right\rvert\right\rvert$$\left(\mathrm{b}\right)$

Instance$\left\lvert\left\lvert\nu_{2}\right\rvert-\left\lvert\nu_{3}\right\rvert\right\rvert$$\left(\mathrm{c}\right)$

Figure 10: This figure depicts absolute differences between (a) 1 and
$\left\lvert\nu_{1}\right\rvert$, (b) $\left\lvert\nu_{1}\right\rvert$ and
$\left\lvert\nu_{2}\right\rvert$, and (c) $\left\lvert\nu_{2}\right\rvert$ and
$\left\lvert\nu_{3}\right\rvert$, where $\left\lvert\nu_{1}\right\rvert$,
$\left\lvert\nu_{2}\right\rvert$ and $\left\lvert\nu_{3}\right\rvert$ denote
the numerically computed first, second and third largest absolute eigenvalues,
respectively, of the transfer matrices of Case 1 in Sec. 4.1. The $2\times
10^{4}$ points on each graph are plotted in ascending order of their
y-coordinate values. The black diamonds denote the maxima of the graphs. The
black circle in (b) represents the minimum of the graph. The black squares in
(c) represent the bounds on the region around $10^{-9}$ where no numerically
sampled differences are observed.

##  Appendix F Values of constants

In this appendix we specify the constants $K$, $c_{1}$, and $c_{2}$ in the
bound (9), as well as the alternative value for $c_{2}$, denoted as $c_{3}$.

First, we determine the values of $c_{1}$ and $K$, and the value of $c_{2}$,
that follow from Theorem III.2 of [SRW15]. We perform the Jordan decomposition
of the transfer matrix $E$,

$E=S\bigoplus_{\nu_{i}}(\nu_{i}\tilde{P}_{\nu_{i}}+\tilde{N}_{\nu_{i}})S^{-1},$
(75)

where $S$ is an invertible matrix, $\nu_{i}$ is an eigenvalue,
$\tilde{P}_{\nu_{i}}$ is the projector onto the subspace corresponding to the
eigenvalue $\nu_{i}$, and $\tilde{N}_{\nu_{i}}$ is a matrix with the elements
immediately above the main diagonal equal to $1$, and all other elements equal
to $0$, possessing the index of nilpotency equal to $(\tilde{K}_{\nu_{i}}+1)$
(i.e., $\tilde{N}_{\nu_{i}}^{\tilde{K}_{\nu_{i}}+1}=0$). Let us set

$\displaystyle K:=\max_{|\nu_{i}|=\nu_{\mathrm{gap}}}\tilde{K}_{\nu_{i}}.$

Let us introduce the set $\mathcal{O}$ containing all distinct values of
$\nu_{i}$.

The Jordan decomposition (75) can contain terms corresponding to equal
eigenvalues. For each distinct eigenvalue, let us combine such the terms
together, arriving to the decomposition

$E=S\bigoplus_{\nu_{i}\in\mathcal{O}}(\nu_{i}P_{\nu_{i}}+N_{\nu_{i}})S^{-1},$
(76)

where $P_{\nu_{i}}$ is still a projector, and $N_{\nu_{i}}$ is still a
nilpotent matrix with the index of nilpotency $K_{\nu_{i}}$ being the largest
one among the nilpotency indices $\tilde{K}_{\nu_{i}}$ of the combined
operators $\tilde{N}_{\nu_{i}}$. We define
$\mathcal{O}_{*}:=\\{\nu_{i}\in\mathcal{O}\;|\;K_{\nu_{i}}=K\\}$, the subset
of $\mathcal{O}$.

Then, by Theorem III.2 in [SRW15], the bound (9) is satisfied for $n$ large
enough that

$\|\bigoplus_{\nu_{i}\in\mathcal{O}}(\nu_{i}P_{\nu_{i}}+N_{\nu_{i}})^{n}\|=\|\max_{\nu_{*}\in\mathcal{O}_{*}}(\nu_{*}P_{\nu_{*}}+N_{\nu_{*}})^{n}\|,$

which may be guaranteed by satisfying the conditions,

$\begin{cases}n>2K,\\\
\left(\frac{n}{K}\right)^{K}\nu_{\mathrm{gap}}^{n-K}\geq(K_{\nu_{i}}+1)\left(\frac{ne}{K_{\nu_{i}}}\right)^{K_{\nu_{i}}}|\nu_{i}|^{n-K_{\nu_{i}}},\mathrm{\;for\;all\;}\nu_{i}\in\mathcal{O}\backslash\mathcal{O}_{*},\end{cases}$

with the constants

$\displaystyle c_{1}$ $\displaystyle:=\|S\|^{-1}\|S^{-1}\|^{-1},$
$\displaystyle c_{2}$
$\displaystyle:=(K+1)\left(\frac{e}{K}\right)^{K}\|S\|\|S^{-1}\|.$

The Theorem IV.3 of [SRW15] provides an alternative bound for
$\|E^{n}-\tilde{E}^{n}\|$,

$\|E^{n}-\tilde{E}^{n}\|\leq\nu_{\mathrm{gap}}^{n+1}\frac{4e^{2}\sqrt{D}(D+1)}{n\left(1-(1+\frac{1}{n})\nu_{\mathrm{gap}}\right)^{\frac{3}{2}}}\sup_{|z|=(1+\frac{1}{n})\nu_{\mathrm{gap}}}\prod_{\nu_{i}\in\mathcal{O}}\left|\frac{1-\bar{\nu_{i}}z}{z-\nu_{i}}\right|^{K_{\nu_{i}}+1},$
(77)

which is satisfied for $n>\nu_{\mathrm{gap}}(1-\nu_{\mathrm{gap}})^{-1}$. Here
$D=\sum_{\nu_{i}\in\mathcal{O}}(K_{\nu_{i}}+1)$, and may be bounded as $D\leq
d_{M}^{2}$. By defining the function

$g(n)=\frac{4e^{2}\sqrt{D}(D+1)}{\left(1-(1+\frac{1}{n})\nu_{\mathrm{gap}}\right)^{\frac{3}{2}}}\left(\frac{\nu_{\mathrm{gap}}}{n}\right)^{K+1}\sup_{|z|=(1+\frac{1}{n})\nu_{\mathrm{gap}}}\prod_{\nu_{i}\in\mathcal{O}}\left|\frac{1-\bar{\nu_{i}}z}{z-\nu_{i}}\right|^{K_{\nu_{i}}+1},$

we may express (77) as

$\|E^{n}-\tilde{E}^{n}\|\leq g(n)n^{K}\nu_{\mathrm{gap}}^{n-K}.$

Notice that $g(n)$ is a bounded function, since for large enough $n$ the
supremum factor in (77) is (i) attained for $z$ in the vanishingly small
neighborhood of some $\nu_{*}\in\mathcal{O}_{*}$, and (ii) grows
asymptotically in proportion to $(n/\nu_{\mathrm{gap}})^{K+1}$. Hence, for
large enough $n$, we may estimate

$\displaystyle\sup_{|z|=(1+\frac{1}{n})\nu_{\mathrm{gap}}}\prod_{\nu_{i}\in\mathcal{O}}\left|\frac{1-\bar{\nu_{i}}z}{z-\nu_{i}}\right|^{K_{\nu_{i}}+1}$
$\displaystyle\leq
2\left(\frac{n}{\nu_{\mathrm{gap}}}\right)^{K+1}\max_{\nu_{*}\in\mathcal{O}_{*}}\prod_{\nu_{*}\neq\nu_{i}\in\mathcal{O}}\left|\frac{1-\bar{\nu_{i}}\nu_{*}}{\nu_{*}-\nu_{i}}\right|^{K_{\nu_{i}}+1},$
(78) $\displaystyle 1-\left(1+\frac{1}{n}\right)\nu_{\mathrm{gap}}$
$\displaystyle>\frac{1-\nu_{\mathrm{gap}}}{2}.$

The precise condition on how large $n$ should be for the first inequality in
(78) to be satisfied is not trivial to specify in a reasonable closed
expression. We choose not to do this, since our primary goal is to simply
illustrate that using the bound of Theorem IV.3 of [SRW15] leads to the
qualitatively same conclusion as using the bound of Theorem III.2 of [SRW15].
With (78), we may bound the $g(n)$ from above,

$\displaystyle
g(n)\leq\frac{16e^{2}\sqrt{D}(D+1)}{\sqrt{2}\left(1-\nu_{\mathrm{gap}}\right)^{\frac{3}{2}}}\left(1-\nu_{\mathrm{gap}}^{2}\right)^{K+1}\left(\frac{2}{\Delta}\right)^{D-K-1},$

where

$\Delta=\min_{\nu_{*}\in\mathcal{O}_{*}}\min_{\nu_{*}\neq\nu_{i}\in\mathcal{O}}|\nu_{i}-\nu_{*}|\leq
2.$

Thus, for large enough $n$, we may bound $g(n)\leq c_{3}$, which leads to
$\|E^{n}-\tilde{E}^{n}\|\leq c_{3}n^{K}\nu_{\mathrm{gap}}^{n-K}$, with

$\displaystyle c_{3}$
$\displaystyle:=\frac{16e^{2}\sqrt{D}(D+1)}{\sqrt{2}\left(1-\nu_{\mathrm{gap}}\right)^{\frac{3}{2}}}\left(1-\nu_{\mathrm{gap}}^{2}\right)^{K+1}\left(\frac{2}{\Delta}\right)^{D-K-1}$
$\displaystyle\leq\frac{16e^{2}d_{M}(d_{M}^{2}+1)}{\sqrt{2}\left(1-\nu_{\mathrm{gap}}\right)^{\frac{3}{2}}}\left(1-\nu_{\mathrm{gap}}^{2}\right)^{K+1}\left(\frac{2}{\Delta}\right)^{d_{M}^{2}-K-1}.$