# $t\geq 1$

Jin Lee Sofía González-García University of California, Santa Barbara
University of California, Santa Barbara<EMAIL_ADDRESS><EMAIL_ADDRESS>Zheng Zhang Haewon Jeong University of
California, Santa Barbara University of California, Santa Barbara
<EMAIL_ADDRESS><EMAIL_ADDRESS>
# Coded Computing Meets Quantum Circuit Simulation: Coded Parallel Tensor
Network Contraction Algorithm

Jin Lee Sofía González-García University of California, Santa Barbara
University of California, Santa Barbara<EMAIL_ADDRESS><EMAIL_ADDRESS>Zheng Zhang Haewon Jeong University of
California, Santa Barbara University of California, Santa Barbara
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Parallel tensor network contraction algorithms have emerged as the pivotal
benchmarks for assessing the classical limits of computation, exemplified by
Google’s demonstration of quantum supremacy through random circuit sampling.
However, the massive parallelization of the algorithm makes it vulnerable to
computer node failures. In this work, we apply coded computing to a practical
parallel tensor network contraction algorithm. To the best of our knowledge,
this is the first attempt to code tensor network contractions. Inspired by
matrix multiplication codes, we provide two coding schemes: 2-node code for
practicality in quantum simulation and hyperedge code for generality. Our
2-node code successfully achieves significant gain for $f$-resilient number
compared to naive replication, proportional to both the number of node
failures and the dimension product of sliced indices. Our hyperedge code can
cover tensor networks out of the scope of quantum, with degraded gain in the
exchange of its generality.

##  I Introduction

The search for quantum advantage hinges on the fact that the dynamics of
quantum computers are hard to simulate with classical resources [1, 2, 3],
given the exponential growth of the state space with the number of qubits.
Tensor networks are arguably the most important classical tool for simulating
quantum circuits. Algorithms such as tensor network contractions [4, 5, 6, 7,
8, 9, 10, 11, 12, 13], together with contraction order optimization [14] and
parallelization by tensor index slicing [12] lie at the heart of determining a
beyond-classical threshold in quantum computer experiments. Google’s quantum
supremacy work relied on this classical benchmark [15, 16]. However, the scope
of tensor network applications extends far beyond quantum circuit simulation.
They are widely used as ansätze for many-body wavefunctions: Matrix Product
States (MPS) for one-dimensional (1d) systems [17, 18], and their higher
dimensional generalization in the form of Projected Entangled Pair States
(PEPS) [19]. Moreover, tensor networks have found potential applications in
machine learning [20, 21, 22, 23], quantum chemistry [24, 25, 26] and other
optimization problems [27, 28].

Despite all the great success and extensive recognition however, tensor
network simulation faces reliability challenges due to the massive
parallelization for efficiency. To effectively simulate state-of-the-art
quantum circuits, tensor network contraction algorithms must be distributed
across millions of classical machines. Distributing workloads to billions of
threads inevitably increases the rate of errors. At the same time, as we are
incorporating more diverse devices with different computational capacities
(e.g., computation-specific accelerators, quantum processors), the variance in
node response time significantly increases. Managing such uncertainties
entirely at the operating systems (OS) level puts too much strain on the
already overloaded OS.

In this paper, we propose a node failure tolerant scheme applied at the
algorithm level by applying coded computing to tensor network contraction.
Coded computing, also known as algorithm-based fault tolerance (ABFT), has
proposed innovative ways to add redundancy in the computation using error-
correcting codes. Coded computing strategies were developed for essential
numerical algorithms such as matrix multiplication [29, 30, 31, 32, 33, 34],
FFT [35, 36], and matrix factorization [37, 38]. Matrix multiplication, in
particular, has been a focal point of extensive research as it is an important
backbone of scientific and machine learning applications.

Our contributions of this paper are for both tensor network and coded
computing community. For the tensor network community, we provide a first work
on parallel tensor network contraction algorithm that can be resilient to
stragglers and failures. For the coded computing community, we provide elegant
generalization of matrix multiplication codes, to the extremely expressive
high-dimensional linear algebraic framework of tensor network.

##  II Background and Notations

In this section, we provide a concise overview of tensor network contraction
and its parallel algorithms, as well as the notations employed throughout the
paper. For a more extensive review of this subject, refer to [39, 40, 41, 42].

###  II-A Tensor Networks and contraction

Tensors are multidimensional arrays that extend the concepts of one-
dimensional vectors and two-dimensional matrices to higher dimensions, and the
tensor network formalism represents linear operations between tensors using
graph structures. A tensor network consists of nodes and edges, where each
node represents a tensor, and each edge connected to the node represents the
index of the tensor. Diagrammatic representations of one-dimensional tensor
$\boldsymbol{A}_{i}$, two-dimensional tensor $\boldsymbol{B}_{i,j}$, and
three-dimensional tensor $\boldsymbol{C}_{i,j,k}$ are shown in Fig 1. We use
bold text with its index labels given in subscript for tensors. E.g.,
$\boldsymbol{B}_{i,j}$ is a two-dimensional tensor with indices labeled as $i$
and $j$; note that this is not a scalar entry $B[i,j]$.
$\boldsymbol{C}_{i,j,k}$ is a three-dimensional tensor and
$\boldsymbol{C}_{1,j,k}$ is a two-dimensional tensor when the first index is
fixed to $i=1$.

An edge connecting two nodes in a tensor network represents a _contraction_
operation, which generalizes an inner product of the dimension—pairwise
multiplication followed by summation. For instance, matrix multiplication can
be represented as Fig. 2(a), where two nodes $\boldsymbol{A}$ and
$\boldsymbol{B}$ are connected through the edge $j$, representing
$C[i,k]=\sum_{j}A[i,j]B[j,k]$. With a slight abuse of notations, we will also
write it in a tensor form as
$\boldsymbol{C}_{i,k}=\sum_{j}\boldsymbol{A}_{i,j}\boldsymbol{B}_{j,k}$. In
Fig. 2(b), we illustrate contracting three edges (indices $j,k,l$) between
three tensors $\boldsymbol{A,B,C}$, resulting in a tensor
$\boldsymbol{T}_{a,b,c}=\sum_{j,k,l}\boldsymbol{A}_{a,j,k}\boldsymbol{B}_{b,k,l}\boldsymbol{C}_{c,l,j}$.
The edges $j,k,l$ are called _closed edges_ and $a,b,c$ are called _open
edges._ Note that the indices of closed edges disappear after the contraction
and the final tensors are indexed only by open edges. With another slight
abuse of notation, we use the terms edge and index interchangeably in this
paper. In tensor networks, there can also be _hyperedges_ that connect more
than two nodes as shown in Fig. 2(c). We call an index shared by $m$ tensors
as an $m$-node index. Based on this diagrammatic representation, we refer an
index shared by more than two tensors as a hyperedge.

When there are multiple closed edges like in Fig. 2(b), contraction is done
sequentially, and the order of contraction does not affect the final outcome.
However, the order of contraction can significantly change the computational
complexity, and choosing an optimal contraction order is a NP-complete
problem[4]. Markov and Shi [4] proved that the computational complexity of the
optimal contraction order for a tensor network is decided by its graph
theoretical property. Hence, the topology of a tensor network decides its
computational cost.

###  II-B Parallel tensor network contraction

In a nutshell, a tensor network is an extremely expressive formalism that
represents a series of arbitrary product sum operations for high-dimensional
tensors. As the tensor network becomes more complex and the sizes of tensors
become large, tensor network contraction can be very computationally
expensive. To scale up the computation, researchers have proposed parallel
tensor network algorithms [12] and showed breakthroughs in quantum circuit
simulation[43, 6, 7].

Figure 1: Diagrammatic representation and corresponding array for
($\boldsymbol{a}$) one dimensional tensor $\boldsymbol{A}_{i}$
($\boldsymbol{b}$) two dimensional tensor $\boldsymbol{B}_{i,j}$
($\boldsymbol{c}$) three dimensional tensor $\boldsymbol{C}_{i,j,k}$ Figure 2:
Diagrammatic representation of tensor network contraction for
($\boldsymbol{a}$) $\sum_{j}\boldsymbol{A}_{i,j}\boldsymbol{B}_{j,k}$
($\boldsymbol{b}$)
$\sum_{k,j,l}\boldsymbol{A}_{a,j,k}\boldsymbol{B}_{b,k,l}\boldsymbol{C}_{c,l,j}$
($\boldsymbol{c}$)
$\sum_{j}\boldsymbol{A}_{j,a}\boldsymbol{B}_{b,j}\boldsymbol{C}_{j,c}$

In this paper, we build up on the work by Chen and Zhang [12], that decomposed
tensor network contraction into an highly parallel algorithm. To explain the
parallelization scheme, we will use an example network given in Figure 3:

$\displaystyle\boldsymbol{T}_{f,g}$
$\displaystyle=\sum_{i,j,k,l,m,n}\boldsymbol{A}_{i}\boldsymbol{B}_{i,j,k}\boldsymbol{C}_{j,l}\boldsymbol{D}_{k,l}\boldsymbol{E}_{k,m}\boldsymbol{F}_{l,n,f}\textbf{G}_{m,n,g}$
$\displaystyle=\sum_{i,j,l,m,n}(\boldsymbol{A}_{i}\boldsymbol{B}_{i,j,1}\boldsymbol{C}_{j,l}\boldsymbol{D}_{1,l}\boldsymbol{E}_{1,m}\boldsymbol{F}_{l,n,f}\textbf{G}_{m,n,g})$
$\displaystyle+\sum_{i,j,l,m,n}(\boldsymbol{A}_{i}\boldsymbol{B}_{i,j,2}\boldsymbol{C}_{j,l}\boldsymbol{D}_{2,l}\boldsymbol{E}_{2,m}\boldsymbol{F}_{l,n,f}\textbf{G}_{m,n,g})$
$\displaystyle\quad\quad\qquad\quad\quad\qquad\quad\quad\qquad\vdots$
$\displaystyle+\sum_{i,j,l,m,n}(\boldsymbol{A}_{i}\boldsymbol{B}_{i,j,L}\boldsymbol{C}_{j,l}\boldsymbol{D}_{L,l}\boldsymbol{E}_{L,m}\boldsymbol{F}_{l,n,f}\textbf{G}_{m,n,g})$
$\displaystyle=\boldsymbol{\sigma}_{1}+\boldsymbol{\sigma}_{2}+\ldots+\boldsymbol{\sigma}_{L}=\boldsymbol{\sigma}_{\text{final}},$
(1)

where $\boldsymbol{\sigma}_{\tilde{k}}$ for $\tilde{k}\in\\{1,2,\cdots,L\\}$
is given as:

$\boldsymbol{\sigma}_{\tilde{k}}=\sum_{i,j,l,m,n}\boldsymbol{A}_{i}\boldsymbol{B}_{i,j,\tilde{k}}\boldsymbol{C}_{j,l}\boldsymbol{D}_{\tilde{k},l}\boldsymbol{E}_{\tilde{k},m}\boldsymbol{F}_{l,n,f}\textbf{G}_{m,n,g}.$
(2)

We say that the tensor network is _sliced_ by the index $k$, and each
$\boldsymbol{\sigma}_{\tilde{k}}$ is _a sliced partition_ of the original
tensor network. Chen and Zhang [12] proposed an algorithm that distributed the
computation of each sliced partition to a separate compute node and summing
them up at the end to obtain $\boldsymbol{\sigma}_{\text{final}}$. For
computing the sliced partition $\boldsymbol{\sigma}_{\tilde{k}}$, no inter-
node communication is required, which makes the algorithm embarrassingly
parallel. In Figure 3, we show the diagrammatic representation of slicing: the
corresponding edge for the slicing index $k$ is erased and the connected
tensors $\boldsymbol{B}_{i,j,k}$, $\boldsymbol{D}_{k,l}$,
$\boldsymbol{E}_{k,m}$ are replaced with sliced subtensors (e.g.,
$\boldsymbol{B}_{i,j,1}$, $\boldsymbol{D}_{1,l}$, $\boldsymbol{E}_{1,m}$)).

For the rest of the paper ,we further simplify the notations by abstracting
out all parts not related to the slicing operation as $\boldsymbol{T}_{-}$,
where $`-^{\prime}$ denotes all indices other than the slicing index. E.g.,
(2) can be written as:

$\boldsymbol{\sigma}_{\tilde{k}}=\sum_{-}(\boldsymbol{B}_{\tilde{k},-}\boldsymbol{D}_{\tilde{k},-}\boldsymbol{E}_{\tilde{k},-})\boldsymbol{T}_{-}.$
(3)

Tensors $\boldsymbol{A,C,F,G}$ that are not connected to the edge $k$ are all
abstracted into a big tensor $\boldsymbol{T}_{-}$. Finally, we say two edges
are _adjacent_ if they are connected with the same node and _non-adjacent_ if
they do not share a node.

Figure 3: Parallelization scheme by slicing index $k$ for the sample tensor
network given in II-B.

##  III System Model and Problem Statement

Figure 4: ($\boldsymbol{a}$) System model for distributed parallel tensor
network contraction without coding. Failure on any computational node will
result as a failure of entire system ($\boldsymbol{b}$) System model for
distributed coded parallel tensor network contraction. Even with some node
failures, as long as the number of successful computational nodes is over
$f$-resilient number, desired final outcome will always be retrieved.

The goal of this paper is to develop coded computing strategies for parallel
tensor network contraction for a given tensor network $G(V,E)$. For parallel
tensor network contraction, we adopt Chen & Zhang’s algorithm [12] described
in Section II-B.We assume a computing model similar to [29, 30] which consists
of a master compute node (_master_) that distributes computation inputs and
assign tasks to worker compute nodes (_workers_). Workers then perform the
assigned computation job and report the result back to the master. With Chen &
Zhang’s algorithm[12], each worker will compute a sliced tensor network
$\boldsymbol{\sigma}_{k}$ and the master node will aggregate the outputs to
obtain $\boldsymbol{\sigma}_{\text{final}}=\sum_{k}\boldsymbol{\sigma}_{k}$ as
shown in Fig. 4.

An important distinction from the model given in [29] is that for high-
dimensional tensor contraction, it is often compute-bound rather than storage-
bound. Thus, instead of a storage constraint, we impose that each worker in
the coded computing scheme cannot have computational costs exceeding the
original uncoded algorithm, i.e., when a worker performs coded computation on
encoded tensor network $\tilde{\boldsymbol{\sigma}}_{k}$, it should satisfy
the following:

$\text{comp. complexity}(\tilde{\boldsymbol{\sigma}}_{k})\leq\text{comp.
complexity}(\boldsymbol{\sigma}_{k}).$ (4)

Computation and memory requirements in tensor network contraction is fully
determined by the network topology111This is assuming that the algorithm does
not optimize for any sparsity in the tensor data., and thus having the same
topology with same tensor dimensions implies the same computation cost.

We assume that in the uncoded algorithm, the tensor network $G$ is sliced
multiple times at $n$ edges, each of which is a $m_{i}$-node edge with its
index dimension $L_{i}$ ($i={1,2,\ldots,n}$). This produces
$N=\prod_{i=1}^{n}L_{i}$ sliced tensor networks and each of them will be
distributed to individual workers. After all $N$ workers complete their
computation, the master aggregates and obtains
$\boldsymbol{\sigma}_{\text{final}}$. We further make an assumption that all
sliced edges are _non-adjacent_ to each other.

We assume that there can be $f$ worker failures and we cannot recover any
partial result from the failed worker. We define _$f$ -resilient number_ as
the number of total workers required in the system to retrieve
$\boldsymbol{\sigma}_{\text{final}}$ in case of arbitrary $f$ node failures.
We will compare our coding strategy with a naive replication strategy, which
has a $f$-resilient number of $N(f+1)$, as we need at least $(f+1)$ replicas
of each worker. We now define our problem below.

Problem Statement. For a tensor network contraction parallelized by slicing
$n$ non-adjacent indices, each of which is $m_{i}$-node and
$L_{i}$-dimensional ($i={1,2,\ldots,n}$), our aim is to design a coded
computing strategy that can be resilient to any $f$ failures while maintaining
the same computational complexity per worker and achieve a $f$-resilient
number smaller than $N(f+1)$ where $N=\prod_{i=1}^{n}L_{i}$.

######  Remark 1.

The assumption that all sliced indices are non-adjacent to each other makes
our problem more tractable and has practical justification. In [6], the
authors suggested a step-dependent iterative heuristics to find indices to
slice by choosing an index that reduces the computational complexity the most
(approximately) at each iteration. With this algorithm, indices connected to
tensors already sliced in previous steps are not likely to be chosen.
Nonetheless, developing a coding scheme for tensor networks allowing adjacent
sliced indices is an interesting open question.

##  IV Main Results

Coded computing strategies we propose in this paper are based on polynomials,
similar to polynomial codes[30] and MatDot/PolyDot codes[29] for distributed
matrix multiplication. I.e., inputs are encoded with polynomials, and after a
set of multilinear computations, output will be encoded in a polynomial of a
higher degree. The $f$-resilient number of the code is determined by the
degree of the resulting polynomial. If the resulting output polynomial is of
$d$-th order, then with any $(d+1)$ successful workers, all its coefficients
can be retrieved by an efficient polynomial interpolation algorithm [44, 30].
For $f$ worker failures, the $f$-resilient number is $(d+f+1)$, and the gain
compared to naive replication is

$\Delta=N(f+1)-(d+f+1)=(N-1)(f+1)-d,$ (5)

where $N=\prod_{i=1}^{n}L_{i}$. Hence, the goal of code design is to minimize
the output polynomial degree while ensuring that all terms relevant to the
tensor network contraction are retrievable from the polynomial coefficients.

In the following subsections, we provide two coding schemes and their
$f$-resilient numbers. The first code is for tensor networks commonly used in
quantum simulation, which are the main interest of contemporary tensor network
studies. The second code is for general tensor networks, applicable for any
tensor network as long as sliced indices are non-adjacent to each other, which
is also a practical assumption.

###  IV-A 2-node Code for Quantum Simulation

While $m$-node hyperedges can represent general tensor operations, in many
practical applications, most edges are conventional 2-node edges. For
instance, tensor networks for quantum simulation often follow the Projected
Entangled-Pair States (PEPS) structure[19], which form a 2D-grid graph that
require only 2-node edges. PEPS is also commonly used in condensed matter
physics [45, 46, 47, 48, 49, 50]. We thus first develop a coding scheme
specifically designed for slicing 2-node indices, which we refer to as the
_2-node code_.

######  Theorem 1.

For a tensor network contraction parallelized by slicing at $n$ non-adjacent
indices, each of which is $2$-node and $L_{i}$-dimensional
($i={1,2,\ldots,n}$), $f$-resilient number of following is achievable:

$f+2\prod_{i=1}^{n}L_{i}-1.$ (6)

The gain compared to replication for $f$ failures is:

$\Delta=(\prod_{i=1}^{n}L_{i}-1)(f-1).$ (7)

###### Proof.

We define template polynomials for $i=1,\cdots,n$:

$\displaystyle p_{i}^{(1)}(x)$ $\displaystyle=1+x+x^{2}+\ldots+x^{L_{i}-1},$
$\displaystyle p_{i}^{(2)}(x)$
$\displaystyle=x^{L_{i}-1}+x^{L_{i}-2}+x^{L_{i}-3}+\ldots+1,$

and then encoding polynomials for index $i$ are
$p_{i}^{(1)}(x^{\prod_{k=1}^{i-1}L_{k}})$ with each coefficient substituted
with the sliced tensors. E.g., for $i=2$, let us call the two tensors
connected to this index as $\mathbf{B}^{(1)}$ and $\mathbf{B}^{(2)}$. Then,
the encoded tensor polynomials are:

$\displaystyle\tilde{\boldsymbol{B}}^{(1)}_{-}(x)$
$\displaystyle=\mathbf{B}^{(1)}_{1,-}+\mathbf{B}^{(1)}_{2,-}x^{L_{1}}+\cdots+\mathbf{B}^{(1)}_{L_{2},-}x^{L_{1}(L_{2}-1)},$
$\displaystyle\tilde{\boldsymbol{B}}^{(2)}_{-}(x)$
$\displaystyle=\mathbf{B}^{(2)}_{1,-}x^{L_{1}(L_{2}-1)}+\mathbf{B}^{(2)}_{2,-}x^{L_{1}(L_{2}-2)}+\cdots+\mathbf{B}^{(2)}_{L_{2},-}.$

Note that all of the desired terms of the product
$\tilde{\boldsymbol{B}}^{(1)}_{-}(x)\tilde{\boldsymbol{B}}^{(2)}_{-}(x)$ are
aligned in the coefficient of $x^{L_{1}(L_{2}-1)}$.

Likewise, the product of encoded tensors of index $i$ is a polynomial of
$2(L_{i}-1)\prod_{k=1}^{i-1}L_{k}$ order which all desired terms are confined
in $(L_{i}-1)\prod_{k=1}^{i-1}L_{k}$ order term. Hence, product of all encoded
tensors is a polynomial of $2((\prod_{i=1}^{n}L_{i})-1)$ order which all
desired terms are confined in $((\prod_{i=1}^{n}L_{i})-1)$ order term. With
any successful $(2(\prod_{i=1}^{n}L_{i})-2+1)$ encoded tensor networks each
using different arbitrary constant $x$, coefficient of
$((\prod_{i=1}^{n}L_{i})-1)$ order term can be retrieved, which is the final
outcome of the original tensor network. Please refer to Appendix A for further
details on the proof. ∎

Additionally, PEPS-structure tensor networks typically consist of uniform
dimensions for all closed edges. Hence, we provide a following corollary for
this special case.

######  Corollary 1.

When $m_{i}=2,L_{i}=L>2$ for all $i$, the following $f$-resilient number is
achievable:

$f+2L^{n}-1.$ (8)

The gain compared to replication for $f$ failures is:

$\Delta=(L^{n}-1)(f-1).$ (9)

######  Example 1 (2-node code).

We provide an example for 2-node code given in Theorem 1 when we have two
sliced edges with $m_{1}=2,L_{1}=4$ and $m_{2}=2,L_{2}=3$. The contraction
equation can be written as:

$\boldsymbol{\sigma}_{a,b}=\sum_{-}(\boldsymbol{A}^{(1)}_{a,-}\boldsymbol{A}^{(2)}_{a,-}\boldsymbol{B}^{(1)}_{b,-}\boldsymbol{B}^{(2)}_{b,-})\boldsymbol{T}_{-},$
(10)

where $a\in\\{1,2,3,4\\}$ and $b\in\\{1,2,3\\}$ represent sliced indices. We
encode the sliced tensors as follows:

$\displaystyle\tilde{\boldsymbol{A}}^{(1)}_{-}$
$\displaystyle=\boldsymbol{A}^{(1)}_{1,-}+\boldsymbol{A}^{(1)}_{2,-}x+\boldsymbol{A}^{(1)}_{3,-}x^{2}+\boldsymbol{A}^{(1)}_{4,-}x^{3},$
$\displaystyle\tilde{\boldsymbol{A}}^{(2)}_{-}$
$\displaystyle=\boldsymbol{A}^{(2)}_{1,-}x^{3}+\boldsymbol{A}^{(2)}_{2,-}x^{2}+\boldsymbol{A}^{(2)}_{3,-}x+\boldsymbol{A}^{(2)}_{4,-},$
$\displaystyle\tilde{\boldsymbol{B}}^{(1)}_{-}$
$\displaystyle=\boldsymbol{B}^{(1)}_{1,-}+\boldsymbol{B}^{(1)}_{2,-}x^{4}+\boldsymbol{B}^{(1)}_{3,-}(x^{4})^{2},$
$\displaystyle\tilde{\boldsymbol{B}}^{(2)}_{-}$
$\displaystyle=\boldsymbol{B}^{(2)}_{1,-}(x^{4})^{2}+\boldsymbol{B}^{(2)}_{2,-}x^{4}+\boldsymbol{B}^{(2)}_{3,-},$

where $x$ can be substituted with a distinct evaluation point at each worker.
After the local contraction operation, the result at each worker is:

$\displaystyle\tilde{\boldsymbol{\sigma}}$
$\displaystyle=\sum_{-}(\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-})\boldsymbol{T}_{-}$
$\displaystyle=(\boldsymbol{\sigma}_{1,1}+\boldsymbol{\sigma}_{2,1}+\boldsymbol{\sigma}_{3,1}+...+\boldsymbol{\sigma}_{4,3})x^{11}+\ldots$
$\displaystyle=\boldsymbol{\sigma}_{\text{final}}x^{11}+\ldots=\text{poly}(x^{22}).$

Hence, with any $(22+1)$ successful workers with distinct evaluation points,
we can perform polynomial interpolation to retrieve
$\boldsymbol{\sigma}_{\text{final}}$. The gain compared to naive replication
is

$\Delta=12(f+1)-(22+f+1)=11(f-1)\hskip 4.0pt(\text{workers}),$

which is consistent with Theorem 1.

###  IV-B Hyperedge Code for General Case

The tensor network framework is increasingly gaining recognition across
various communities outside of the quantum field, and tensor network
structures that incorporate hyperedges may find more use in the future. We
thus provide a more general coding scheme beyond the 2-node code, suitable for
any tensor network, which we refer to as the _hyperedge code_. We establish an
achievable scheme for hyperedge codes in Theorem 2. The provided code works
for any tensor network as long as two sliced edges are not adjacent to each
other.

######  Theorem 2.

For a tensor network contraction parallelized by slicing at $n$ non-adjacent
indices, each of which is $m_{i}$-node and $L_{i}$-dimensional
($i={1,2,\ldots,n}$), $f$-resilient number of following is achievable:

$f+\prod_{i=1}^{n}\frac{(m^{L_{i}}_{i}-1)}{(m_{i}-1)}.$ (11)

The corresponding gain compared to naive replication is

$\Delta=(\prod_{i=1}^{n}L_{i}-1)(f+1)-\prod_{i=1}^{n}\frac{(m^{L_{i}}_{i}-1)}{(m_{i}-1)}+1.$
(12)

###### Proof.

We define template polynomial for $i=1,\cdots,n$:

$\displaystyle
p_{i}(x)=1+x+x^{1+m_{i}}+x^{1+m_{i}+m_{i}^{2}}+\ldots+x^{1+\dots+m_{i}^{L_{i}-2}},$

and then the encoding polynomial for index $i$ is
$p_{i}(x^{\prod_{j=1}^{i-1}\frac{m_{j}^{L_{j}}-1}{m_{j}-1}})$ with each
coefficient substituted with the sliced tensors. E.g., for $i=2$, let us call
the $m_{2}$ tensors connected to this index as $\mathbf{B}^{(j)}$ where
$j=1,\ldots,m_{2}$. Then, the encoded tensor polynomials are:

$\displaystyle\tilde{\boldsymbol{B}}^{(j)}_{-}(x)=$
$\displaystyle\mathbf{B}^{(j)}_{1,-}+\mathbf{B}^{(j)}_{2,-}x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}+\cdots$
$\displaystyle+\mathbf{B}^{(j)}_{L_{2},-}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})^{1+\ldots+m_{2}^{L_{2}-2}},$

Note that every desired term of the product
$\prod_{j=1}^{m_{2}}\tilde{\boldsymbol{B}}^{(j)}_{-}(x)$ is each confined in
different order term only by itself, and the resulting polynomial is
$\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\cdot\frac{m_{2}^{L_{2}}-m_{2}}{m_{2}-1}$
order. Likewise, the product of encoded tensors of index $i$ is a polynomial
of
$\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\cdot\frac{m_{2}^{L_{2}}-1}{m_{2}-1}\cdot\frac{m_{3}^{L_{3}}-1}{m_{3}-1}\cdot...\cdot\frac{m_{i}^{L_{i}}-m_{i}}{m_{i}-1}$
order which each desired term is confined in different order term only by
itself. Hence, product of all encoded tensors is a polynomial of
$(\prod_{j=1}^{n}\frac{m_{j}^{L_{j}}-1}{m_{j}-1})-1$ order which each desired
term is confined in different order term only by itself. With any successful
$\prod_{j=1}^{n}\frac{m_{j}^{L_{j}}-1}{m_{j}-1}$ encoded tensor networks each
using different arbitrary constant $x$, all the desired outcomes of sliced
tensor networks can be retrieved and the final outcome of the original tensor
network can be acquired by their summation. Please refer to Appendix B for
further details on the proof.∎

######  Example 2 (Hyperedge code).

We provide an example for hyperedge code given in Theorem 2 when we have two
sliced edges with $m_{1}=3,L_{1}=2$ and $m_{2}=4,L_{2}=2$. The contraction
equation can be written as:

$\boldsymbol{\sigma}_{a,b}=\sum_{-}(\boldsymbol{A}^{(1)}_{a,-}\boldsymbol{A}^{(2)}_{a,-}\boldsymbol{A}^{(3)}_{a,-}\boldsymbol{B}^{(1)}_{b,-}\boldsymbol{B}^{(2)}_{b,-}\boldsymbol{B}^{(3)}_{b,-}\boldsymbol{B}^{(4)}_{b,-})\boldsymbol{T}_{-},$

where $a\in\\{1,2\\}$ and $b\in\\{1,2\\}$ represent sliced indices. We encode
the sliced tensors as

$\displaystyle\tilde{\boldsymbol{A}}^{(i)}_{-}$
$\displaystyle=\boldsymbol{A}^{(i)}_{1,-}+\boldsymbol{A}^{(i)}_{2,-}x$
$\displaystyle\tilde{\boldsymbol{B}}^{(j)}_{-}$
$\displaystyle=\boldsymbol{B}^{(j)}_{1,-}+\boldsymbol{B}^{(j)}_{2,-}x^{1+3}$

where $i=1,2,3$ and $j=1,2,3,4$. $x$ can be substituted with a distinct
evaluation point at each worker. After the local contraction operation, the
result at each worker is:

$\displaystyle\tilde{\boldsymbol{\sigma}}$
$\displaystyle=\sum_{-}(\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}\tilde{\boldsymbol{A}}^{(3)}_{-}\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-}\tilde{\boldsymbol{B}}^{(3)}_{-}\tilde{\boldsymbol{B}}^{(4)}_{-})\boldsymbol{T}_{-}$
$\displaystyle=\boldsymbol{\sigma}_{1,1}+\boldsymbol{\sigma}_{2,1}x^{3}+\boldsymbol{\sigma}_{1,2}x^{4}+\boldsymbol{\sigma}_{2,2}x^{19}+\ldots$
$\displaystyle=\text{poly}(x^{19}).$

Hence, with any $(19+1)$ successful workers with distinct evaluation points,
we can perform polynomial interpolation to retrieve
$\boldsymbol{\sigma}_{\text{final}}$. The gain compared to naive replication
is

$\Delta=4(f+1)-(f+20)=3f-16\hskip 4.0pt(\text{workers}),$

which is consistent with Theorem 2.

######  Remark 2.

Negative product term of $m_{i}$ in (12) of Theorem 2 degrades the gain of the
hyperedge code. Intuitively, that is because hyperedge code confines all
desired terms each in different order terms hence resulting polynomial
requires higher degree for larger $\\{m_{i}\\}$. On the contrary, 2-node code
confines all desired terms in a single term of the resulting polynomial,
therefore its gain (7) does not suffer such loss. However, we emphasize that
our codes and their upper bounds of $f$-resilient number are not proven to be
optimal nor tight.

###  IV-C Master Node Complexity

Lastly, we show that encoding and decoding costs are negligible compared to
computation-intensive tensor contraction. We will make an informal argument
with as simple setting where all edges are $m$-node and $L$-dimensional, and
all nodes have a degree $k$. The computation cost of tensor contraction is
determined by the most expensive edge contraction [4]. In our simplified
setting, all edge contraction has the same cost. Contracting any index will
result in an intermediate tensor with $m(k-1)$ open edges, as there are $m$
tensors with each $(k-1)$ edges not being contracted. The computational
complexity for this is $O(L^{(k-1)m}\cdot L)$. On the other hand, encoding
sliced tensors is equivalent to taking a linear combination of $L$ sliced
tensors each of which has $k-1$ open edges. Hence, the cost of encoding is
$O(L^{k-1}\cdot L)$, assuming $m$ is not exponentially large compared to $L$.
Compared to contraction complexity of $O(L^{(k-1)m+1})$, encoding complexity
of $O(L^{k})$ is negligible. For example, for a typical PEPS tensor network
with $k=4$ and $m=2$, contraction cost for an edge is $O(L^{7})$ and encoding
cost is $O(L^{4})$. For decoding complexity, the final outcome of tensor
network contraction in most quantum applications is a scalar or a vector, and
thus decoding is minimal [44, 30].

###  IV-D Existing matrix codes under our framework

Figure 5: Diagrammatic representation of ($\boldsymbol{a}$) polynomial code
and ($\boldsymbol{b}$) MatDot code in tensor network formalism.

As the tensor network formalism is a general linear algebraic framework that
also encompasses matrix multiplication, our work can be thought of as a
generalized expansion of precedent matrix multiplication coding schemes. For
instance, in the tensor network formalism, MatDot codes[29] slice a 2-node
index between $\boldsymbol{A}$ and $\boldsymbol{B}$ (see Fig. 5), which
corresponds to the case of $n=1$, $m_{1}=2$, and $N=L_{1}$ for Theorem 1. In
the formalism of our work, $f$-resilient number of the 2-node code is
$f+2L_{1}-1$, which is consistent with the $f$-resilient number of MatDot
codes [29]. On the other hand, polynomial code[30] corresponds to slicing
through open edges in tensor network formalism (see Fig. 5). This is a case we
did not discuss in our work as open edges are not to be sliced as they do not
reduce the complexity of the tensor network. Finally, $n$-matrix
multiplication codes described in [29] correspond to a chain in tensor network
formalism and slicing through all edges (with alternating index dimension of
$s$ and $t$). While this gave a rise to high $f$-resilient number in the
original work[29], our work shines a new light on this problem that in
practical algorithms, coding overhead can be much lower than $n$-matrix codes
as we only slice through a select few edges instead of all edges [12].

##  V Discussion and future works

This work is the first attempt at applying coded computing for reliable tensor
network contraction on large-scale parallel machines. There are many more
exciting future directions to explore. While we provide a coding scheme that
works only for non-adjacent edges with some practical justification, making it
more general to include both adjacent and non-adjacent edge slicing would be
an interesting problem. Another crucial research question is a lower bound on
the $f$-resilient number of coded tensor network algorithms. While we provide
two different achievable schemes in Theorem 1 and Theorem 2, we do not know if
the proposed approach is optimal or how far from optimality it is. Finally, as
we are inspired by very practical applications, it would be valuable to
implement the proposed coded computing strategies and perform an empirical
evaluation on real-world high-performance computing systems.

## References

  * [1] M. J. Bremner, R. Jozsa, and D. J. Shepherd, “Classical simulation of commuting quantum computations implies collapse of the polynomial hierarchy,” _Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences_ , vol. 467, no. 2126, p. 459–472, Aug. 2010. [Online]. Available: http://dx.doi.org/10.1098/rspa.2010.0301
  * [2] S. Aaronson and A. Arkhipov, “The computational complexity of linear optics,” 2010.
  * [3] M. J. Bremner, A. Montanaro, and D. J. Shepherd, “Average-case complexity versus approximate simulation of commuting quantum computations,” _Physical Review Letters_ , vol. 117, no. 8, Aug. 2016. [Online]. Available: http://dx.doi.org/10.1103/PhysRevLett.117.080501
  * [4] I. L. Markov and Y. Shi, “Simulating quantum computation by contracting tensor networks,” _SIAM Journal on Computing_ , vol. 38, no. 3, pp. 963–981, 2008.
  * [5] S. Boixo, S. V. Isakov, V. N. Smelyanskiy, and H. Neven, “Simulation of low-depth quantum circuits as complex undirected graphical models,” _arXiv preprint arXiv:1712.05384_ , 2017.
  * [6] D. Lykov, R. Schutski, A. Galda, V. Vinokur, and Y. Alexeev, “Tensor network quantum simulator with step-dependent parallelization,” in _2022 IEEE International Conference on Quantum Computing and Engineering (QCE)_ , 2022, pp. 582–593.
  * [7] T. Vincent, L. J. O’Riordan, M. Andrenkov, J. Brown, N. Killoran, H. Qi, and I. Dhand, “Jet: Fast quantum circuit simulations with parallel task-based tensor-network contraction,” _Quantum_ , vol. 6, p. 709, 2022.
  * [8] Y. Liu, Y. Chen _et al._ , “Verifying quantum advantage experiments with multiple amplitude tensor network contraction,” _Physical Review Letters_ , vol. 132, no. 3, Jan. 2024. [Online]. Available: http://dx.doi.org/10.1103/PhysRevLett.132.030601
  * [9] C. Huang, F. Zhang, M. Newman, J. Cai, X. Gao, Z. Tian, J. Wu, H. Xu, H. Yu, B. Yuan, M. Szegedy, Y. Shi, and J. Chen, “Classical simulation of quantum supremacy circuits,” 2020.
  * [10] F. Pan, K. Chen, and P. Zhang, “Solving the sampling problem of the sycamore quantum circuits,” _Physical Review Letters_ , vol. 129, no. 9, 2022.
  * [11] G. Kalachev, P. Panteleev, P. Zhou, and M.-H. Yung, “Classical sampling of random quantum circuits with bounded fidelity,” 2021.
  * [12] J. Chen, F. Zhang, C. Huang, M. Newman, and Y. Shi, “Classical simulation of intermediate-size quantum circuits,” _arXiv preprint arXiv:1805.01450_ , 2018.
  * [13] B. Villalonga, S. Boixo, B. Nelson, C. Henze, E. Rieffel, R. Biswas, and S. Mandrà, “A flexible high-performance simulator for verifying and benchmarking quantum circuits implemented on real hardware,” _npj Quantum Information_ , vol. 5, no. 1, Oct. 2019. [Online]. Available: http://dx.doi.org/10.1038/s41534-019-0196-1
  * [14] J. Gray and S. Kourtis, “Hyper-optimized tensor network contraction,” _Quantum_ , vol. 5, p. 410, Mar. 2021. [Online]. Available: http://dx.doi.org/10.22331/q-2021-03-15-410
  * [15] F. Arute, K. Arya, R. Babbush, D. Bacon, J. C. Bardin, R. Barends, R. Biswas, S. Boixo, F. G. Brandao, D. A. Buell _et al._ , “Quantum supremacy using a programmable superconducting processor,” _Nature_ , vol. 574, no. 7779, pp. 505–510, 2019.
  * [16] A. Morvan, B. Villalonga, X. Mi, S. Mandrà _et al._ , “Phase transition in random circuit sampling,” 2023.
  * [17] S. R. White, “Density matrix formulation for quantum renormalization groups,” _Phys. Rev. Lett._ , vol. 69, pp. 2863–2866, Nov 1992. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevLett.69.2863
  * [18] M. Fannes, B. Nachtergaele, and R. F. Werner, “Finitely correlated states on quantum spin chains,” _Communications in mathematical physics_ , vol. 144, pp. 443–490, 1992.
  * [19] F. Verstraete and J. I. Cirac, “Renormalization algorithms for quantum-many body systems in two and higher dimensions,” 2004.
  * [20] E. M. Stoudenmire and D. J. Schwab, “Supervised learning with quantum-inspired tensor networks,” 2017.
  * [21] W. Huggins, P. Patil, B. Mitchell, K. B. Whaley, and E. M. Stoudenmire, “Towards quantum machine learning with tensor networks,” _Quantum Science and Technology_ , vol. 4, no. 2, p. 024001, Jan. 2019. [Online]. Available: http://dx.doi.org/10.1088/2058-9565/aaea94
  * [22] I. Glasser, N. Pancotti, and J. I. Cirac, “From probabilistic graphical models to generalized tensor networks for supervised learning,” 2019.
  * [23] S. Cheng, L. Wang, and P. Zhang, “Supervised learning with projected entangled pair states,” _Physical Review B_ , vol. 103, no. 12, p. 125117, 2021.
  * [24] D.-T. Chen, P. Helms, A. R. Hale, M. Lee, C. Li, J. Gray, G. Christou, V. S. Zapf, G. K.-L. Chan, and H.-P. Cheng, “Using hyper-optimized tensor networks and first-principles electronic structure to simulate experimental properties of the giant Mn84 torus,” 2022.
  * [25] K. H. Marti and M. Reiher, “The density matrix renormalization group algorithm in quantum chemistry,” _Zeitschrift für Physikalische Chemie_ , vol. 224, no. 3-4, pp. 583–599, 2010. [Online]. Available: https://doi.org/10.1524/zpch.2010.6125
  * [26] S. Lee, J. Lee, H. Zhai, Y. Tong, A. M. Dalzell, A. Kumar, P. Helms, J. Gray, Z.-H. Cui, W. Liu, M. Kastoryano, R. Babbush, J. Preskill, D. R. Reichman, E. T. Campbell, E. F. Valeev, L. Lin, and G. K.-L. Chan, “Evaluating the evidence for exponential quantum advantage in ground-state quantum chemistry,” _Nature Communications_ , vol. 14, no. 1, Apr. 2023. [Online]. Available: http://dx.doi.org/10.1038/s41467-023-37587-6
  * [27] S. Kourtis, C. Chamon, E. Mucciolo, and A. Ruckenstein, “Fast counting with tensor networks,” _SciPost Physics_ , vol. 7, no. 5, Nov. 2019. [Online]. Available: http://dx.doi.org/10.21468/SciPostPhys.7.5.060
  * [28] J.-G. Liu, X. Gao, M. Cain, M. D. Lukin, and S.-T. Wang, “Computing solution space properties of combinatorial optimization problems via generic tensor networks,” _SIAM Journal on Scientific Computing_ , vol. 45, no. 3, p. A1239–A1270, Jun. 2023. [Online]. Available: http://dx.doi.org/10.1137/22M1501787
  * [29] S. Dutta, M. Fahim, F. Haddadpour, H. Jeong, V. Cadambe, and P. Grover, “On the optimal recovery threshold of coded matrix multiplication,” _IEEE Transactions on Information Theory_ , vol. 66, no. 1, pp. 278–301, 2020.
  * [30] Q. Yu, M. Maddah-Ali, and S. Avestimehr, “Polynomial codes: an optimal design for high-dimensional coded matrix multiplication,” _Advances in Neural Information Processing Systems_ , vol. 30, 2017.
  * [31] S. Wang, J. Liu, and N. Shroff, “Coded sparse matrix multiplication,” in _International Conference on Machine Learning (ICML)_ , 2018, pp. 5139–5147.
  * [32] A. Severinson, A. G. i Amat, and E. Rosnes, “Block-Diagonal and LT Codes for Distributed Computing With Straggling Servers,” _IEEE Transactions on Communications_ , vol. 67, no. 3, pp. 1739–1753, 2019.
  * [33] A. Mallick, M. Chaudhari, and G. Joshi, “Fast and efficient distributed matrix-vector multiplication using rateless fountain codes,” in _IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP)_ , 2019, pp. 8192–8196.
  * [34] H. Jeong, A. Devulapalli, V. R. Cadambe, and F. P. Calmon, “e-approximate coded matrix multiplication is nearly twice as efficient as exact multiplication,” _IEEE Journal on Selected Areas in Information Theory_ , vol. 2, no. 3, pp. 845–854, 2021.
  * [35] H. Jeong, T. M. Low, and P. Grover, “Masterless Coded Computing: A Fully-Distributed Coded FFT Algorithm,” in _IEEE Communication, Control, and Computing (Allerton)_ , 2018, pp. 887–894.
  * [36] Q. Yu, M. A. Maddah-Ali, and A. S. Avestimehr, “Coded fourier transform,” in _2017 55th Annual Allerton Conference on Communication, Control, and Computing (Allerton)_. IEEE, 2017, pp. 494–501.
  * [37] Q. M. Nguyen, H. Jeong, and P. Grover, “Coded QR Decomposition,” in _IEEE International Symposium on Information Theory (ISIT)_ , 2020.
  * [38] Q. M. Nguyen, I. Weissburg, and H. Jeong, “Coded computing for fault-tolerant parallel qr decomposition,” 2023.
  * [39] G. Evenbly, “A practical guide to the numerical implementation of tensor networks i: Contractions, decompositions, and gauge freedom,” _Frontiers in Applied Mathematics and Statistics_ , vol. 8, p. 806549, 2022.
  * [40] J. Biamonte and V. Bergholm, “Tensor networks in a nutshell,” _arXiv preprint arXiv:1708.00006_ , 2017.
  * [41] R. Orús, “Tensor networks for complex quantum systems,” _Nature Reviews Physics_ , vol. 1, no. 9, pp. 538–550, 2019.
  * [42] S. Montangero, E. Montangero, and Evenson, _Introduction to tensor network methods_. Springer, 2018.
  * [43] C. Huang, F. Zhang, M. Newman, X. Ni, D. Ding, J. Cai, X. Gao, T. Wang, F. Wu, G. Zhang _et al._ , “Efficient parallelization of tensor network contraction for simulating quantum computation,” _Nature Computational Science_ , vol. 1, no. 9, pp. 578–587, 2021.
  * [44] K. S. Kedlaya and C. Umans, “Fast polynomial factorization and modular composition,” _SIAM Journal on Computing_ , vol. 40, no. 6, pp. 1767–1802, 2011.
  * [45] P. Corboz and F. Mila, “Crystals of bound states in the magnetization plateaus of the shastry-sutherland model,” _Physical Review Letters_ , vol. 112, no. 14, Apr. 2014. [Online]. Available: http://dx.doi.org/10.1103/PhysRevLett.112.147203
  * [46] P. Corboz, T. M. Rice, and M. Troyer, “Competing states in the t-j model: Uniform d-wave state versus stripe state,” _Physical review letters_ , vol. 113, no. 4, p. 046402, 2014.
  * [47] I. Niesen and P. Corboz, “A tensor network study of the complete ground state phase diagram of the spin-1 bilinear-biquadratic Heisenberg model on the square lattice,” _SciPost Phys._ , vol. 3, p. 030, 2017. [Online]. Available: https://scipost.org/10.21468/SciPostPhys.3.4.030
  * [48] B.-X. Zheng, C.-M. Chung, P. Corboz, G. Ehlers, M.-P. Qin, R. M. Noack, H. Shi, S. R. White, S. Zhang, and G. K.-L. Chan, “Stripe order in the underdoped region of the two-dimensional hubbard model,” _Science_ , vol. 358, no. 6367, p. 1155–1160, Dec. 2017. [Online]. Available: http://dx.doi.org/10.1126/science.aam7127
  * [49] B. Ponsioen, S. S. Chung, and P. Corboz, “Period 4 stripe in the extended two-dimensional hubbard model,” _Physical Review B_ , vol. 100, no. 19, Nov. 2019. [Online]. Available: http://dx.doi.org/10.1103/PhysRevB.100.195141
  * [50] J.-Y. Chen, S. Capponi, A. Wietek, M. Mambrini, N. Schuch, and D. Poilblanc, “Su (3) 1 chiral spin liquid on the square lattice: A view from symmetric projected entangled pair states,” _Physical Review Letters_ , vol. 125, no. 1, p. 017201, 2020.

##  Appendix A Proof of Theorem 1

Let us refer to the sliced indices by number $i={1,2,\ldots,n}$. Each index
$i$ is $L_{i}$ dimensional and shared by 2 tensors. For convenience, let us
denote the tensors sharing index $1$ as
$(\boldsymbol{A}^{(1)}_{-},\boldsymbol{A}^{(2)}_{-})$, the tensors sharing
index $2$ as $(\boldsymbol{B}^{(1)}_{-},\boldsymbol{B}^{(2)}_{-})$, an the
tensors sharing index $n$ as
$(\boldsymbol{Z}^{(1)}_{-},\boldsymbol{Z}^{(2)}_{-})$. Each of the tensors is
sliced by the selected index, for instance, $\boldsymbol{A}^{(1)}_{-}$ is
sliced as $\boldsymbol{A}^{(1)}_{s_{1},-}$ for fixed value
$s_{1}\in\\{1,2,\ldots,L_{1}\\}$ of index $1$. Then for fixed value $s_{i}$
for each index $i$, which span by $s_{i}\in\\{1,2,\ldots,L_{i}\\}$, equation
for each sliced tensor network corresponding to the string
$s_{1}s_{2}...s_{n}$ is

$\boldsymbol{\sigma}_{s_{1}s_{2}...s_{n}}=\sum_{-}(\boldsymbol{A}^{(1)}_{s_{1},-}\boldsymbol{A}^{(2)}_{s_{1},-}\boldsymbol{B}^{(1)}_{s_{2},-}\boldsymbol{B}^{(2)}_{s_{2},-}...\boldsymbol{Z}^{(1)}_{s_{n},-}\boldsymbol{Z}^{(2)}_{s_{n},-})\boldsymbol{T}_{-},$

where each index $i$ spans by $s_{i}\in\\{1,2,\ldots,L_{i}\\}$. Hence, there
are $(N=\prod_{i=1}^{n}L_{i})$ terms in the summation we need to retrieve.

For index $i$, we define two template polynomials:

$p_{i}^{(1)}(x)=1+x+x^{2}+\ldots+x^{L_{i}-1}$
$p_{i}^{(2)}(x)=x^{L_{i}-1}+x^{L_{i}-2}+x^{L_{i}-3}+\ldots+1.$

We encode for index $1$ with encoding polynomial $p_{1}^{(1)}(x)$ and
$p_{2}^{(2)}(x)$ by putting each sliced tensors of $\boldsymbol{A}^{(1)}_{-}$
and $\boldsymbol{A}^{(2)}_{-}$ to the coefficients of $p_{1}^{(1)}(x)$ and
$p_{2}^{(2)}(x)$ respectively as below:

$\tilde{\boldsymbol{A}}^{(1)}_{-}=\boldsymbol{A}^{(1)}_{1,-}+x\boldsymbol{A}^{(1)}_{2,-}+\ldots+x^{L_{1}-1}\boldsymbol{A}^{(1)}_{L_{1},-}=\sum_{k=1}^{L_{1}}x^{k-1}\boldsymbol{A}^{(1)}_{k,-}$
$\tilde{\boldsymbol{A}}^{(2)}_{-}=x^{L_{1}-1}\boldsymbol{A}^{(2)}_{1,-}+x^{L_{1}-2}\boldsymbol{A}^{(2)}_{2,-}+..+\boldsymbol{A}^{(2)}_{L_{1},-}=\sum_{k=1}^{L_{1}}x^{L_{1}-k}\boldsymbol{A}^{(2)}_{k,-}.$

Then, product of encoded tensors
$\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}$ is:

$\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}=(\sum_{k=1}^{L_{1}}x^{k-1}\boldsymbol{A}^{(1)}_{k,-})(\sum_{k=1}^{L_{1}}x^{L_{1}-k}\boldsymbol{A}^{(2)}_{k,-})$
$=x^{L_{1}-1}(\boldsymbol{A}^{(1)}_{1,-}\boldsymbol{A}^{(2)}_{1,-}+\boldsymbol{A}^{(1)}_{2,-}\boldsymbol{A}^{(2)}_{2,-}+\ldots+\boldsymbol{A}^{(1)}_{L_{1},-}\boldsymbol{A}^{(2)}_{L_{1},-})$
$+(\text{other order terms})=poly(x^{2(L_{1}-1)}).$

All the desired terms of
$\boldsymbol{A}^{(1)}_{s_{1},-}\boldsymbol{A}^{(2)}_{s_{1},-}(s_{1}\in\\{1,2,...,L_{1}\\})$
are confined in one coefficient of $(L_{1}-1)$ order term.

For index $2$, we use encoding polynomial $p^{(1)}_{2}(x^{L_{1}})$ and
$p^{(2)}_{2}(x^{L_{1}})$:

$\tilde{\boldsymbol{B}}^{(1)}_{-}=\boldsymbol{B}^{(1)}_{1,-}+\ldots+x^{L_{1}(L_{2}-1)}\boldsymbol{B}^{(1)}_{L_{2},-}=\sum_{k=1}^{L_{2}}x^{L_{1}(k-1)}\boldsymbol{B}^{(1)}_{k,-}$
$\tilde{\boldsymbol{B}}^{(2)}_{-}=x^{L_{1}(L_{2}-1)}\boldsymbol{B}^{(2)}_{1,-}+\ldots+\boldsymbol{B}^{(2)}_{L_{2},-}=\sum_{k=1}^{L_{2}}x^{L_{1}(L_{2}-k)}\boldsymbol{B}^{(2)}_{k,-}.$
$\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-}=(\sum_{k=1}^{L_{2}}x^{L_{1}(k-1)}\boldsymbol{B}^{(1)}_{k,-})(\sum_{k=1}^{L_{2}}x^{L_{1}(L_{2}-k)}\boldsymbol{B}^{(2)}_{k,-})$
$=x^{L_{1}(L_{2}-1)}(\boldsymbol{B}^{(1)}_{1,-}\boldsymbol{B}^{(2)}_{1,-}+\boldsymbol{B}^{(1)}_{2,-}\boldsymbol{B}^{(2)}_{2,-}+\ldots+\boldsymbol{B}^{(1)}_{L_{2},-}\boldsymbol{B}^{(2)}_{L_{2},-})$
$+(\text{other order terms})=poly(x^{2L_{1}(L_{2}-1)}).$

Likewise, all the desired terms of
$\boldsymbol{B}^{(1)}_{s_{2},-}\boldsymbol{B}^{(2)}_{s_{2},-}(s_{2}\in\\{1,2,...,L_{2}\\})$
are confined in one coefficient of $L_{1}(L_{2}-1)$ order term.

Furthermore, polynomial resulting from the product
($\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-}$)
is:

$\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-}=x^{L_{1}L_{2}-1}(\sum_{s_{1},s_{2}}\boldsymbol{A}^{(1)}_{s_{1},-}\boldsymbol{A}^{(2)}_{s_{1},-}\boldsymbol{B}^{(1)}_{s_{2}-}\boldsymbol{B}^{(2)}_{s_{2},-})$
$+(\text{other order terms})=poly(x^{2(L_{1}L_{2}-1)}).$

All desired terms of
$(\boldsymbol{A}^{(1)}_{s_{1},-}\boldsymbol{A}^{(2)}_{s_{1},-}\boldsymbol{B}^{(1)}_{s_{2}-}\boldsymbol{B}^{(2)}_{s_{2},-})$
for $s_{1}\in\\{1,2,..,L_{1}\\}$ $s_{2}\in\\{1,2,..,L_{2}\\}$ are confined in
one coefficient of $(L_{1}L_{2}-1)$ order term.

For index $3$, encoding polynomials are $p^{(1)}_{3}(x^{L_{1}L_{2}})$ and
$p^{(2)}_{3}(x^{L_{1}L_{2}})$. For index $i$, encoding polynomials are
$p^{(1)}_{i}(x^{\prod_{j=1}^{i-1}L_{j}})$ and
$p^{(2)}_{i}(x^{\prod_{j=1}^{i-1}L_{j}})$. For the last index $n$, encoded
tensors are

$\tilde{\boldsymbol{Z}}^{(1)}_{-}=\sum_{k=1}^{L_{n}}x^{(k-1)\prod_{j=1}^{n-1}L_{j}}\boldsymbol{Z}^{(1)}_{k,-}$
$\tilde{\boldsymbol{Z}}^{(2)}_{-}=\sum_{k=1}^{L_{n}}x^{(L_{n}-k)\prod_{j=1}^{n-1}L_{j}}\boldsymbol{Z}^{(2)}_{k,-}.$

And their product is

$\tilde{\boldsymbol{Z}}^{(1)}_{-}\tilde{\boldsymbol{Z}}^{(2)}_{-}=x^{(L_{n}-1)\prod_{j=1}^{n-1}L_{j}}(\boldsymbol{Z}^{(1)}_{1,-}\boldsymbol{Z}^{(2)}_{1,-}+\ldots+\boldsymbol{Z}^{(1)}_{L_{n},-}\boldsymbol{Z}^{(2)}_{L_{n},-}).$

As a result, an encoded tensor network produces a polynomial at the end of the
contraction, which confines all desired outcome of sliced tensor networks in
one coefficient of $((\prod_{i=1}^{n}L_{i})-1)$ order term.

$\tilde{\boldsymbol{\sigma}}=\sum_{-}(\tilde{\boldsymbol{A}}^{(1)}_{-}\tilde{\boldsymbol{A}}^{(2)}_{-}\tilde{\boldsymbol{B}}^{(1)}_{-}\tilde{\boldsymbol{B}}^{(2)}_{-}...\tilde{\boldsymbol{Z}}^{(1)}_{-}\tilde{\boldsymbol{Z}}^{(2)}_{-})\boldsymbol{T}_{-}$
$=\sum_{-}[(\sum_{s_{1},s_{2},..,s_{n}}\boldsymbol{A}^{(1)}_{s_{1},-}\boldsymbol{A}^{(2)}_{s_{1},-}\boldsymbol{B}^{(1)}_{s_{2},-}..\boldsymbol{Z}^{(2)}_{s_{n},-})x^{(\prod_{i=1}^{n}L_{i})-1}+..]\boldsymbol{T}_{-}$
$=(\sum_{s_{1},s_{2},..,s_{n}}\boldsymbol{\sigma}_{s_{1}s_{2}...s_{n}})x^{(\prod_{i=1}^{n}L_{i})-1}+(\text{other
order terms})$
$=\boldsymbol{\sigma}_{\text{final}}x^{(\prod_{i=1}^{n}L_{i})-1}+(\text{other
order terms})=poly(x^{2(\prod_{i=1}^{n}L_{i})-2}).$

Hence, with any successful $(2(\prod_{i=1}^{n}L_{i})-2+1)$ encoded tensor
networks each using different arbitrary constant $x$, coefficient of
$((\prod_{i=1}^{n}L_{i})-1)$ order term in the resulting polynomial can be
retrieved, which is the final outcome of the original tensor network.

Therefore, for $f$ node failures, $f$-resilient number is

$f+2\prod_{i=1}^{n}L_{i}-1$

for 2-node code.

##  Appendix B Proof of Theorem 2

Let us refer to the sliced indices by number $i={1,2,\ldots,n}$. Each index
$i$ is $L_{i}$ dimensional and shared by $m_{i}$ tensors. For convenience, let
us denote the tensors sharing index $1$ as
$(\boldsymbol{A}^{(1)}_{-},\boldsymbol{A}^{(2)}_{-},\ldots,\boldsymbol{A}^{(m_{1})}_{-})$,
and tensors sharing index $2$ as
$(\boldsymbol{B}^{(1)}_{-},\boldsymbol{B}^{(2)}_{-},\ldots,\boldsymbol{B}^{(m_{2})}_{-})$.
We will denote tensors sharing index $n$ as
$(\boldsymbol{Z}^{(1)}_{-},\boldsymbol{Z}^{(2)}_{-},\ldots,\boldsymbol{Z}^{(m_{n})}_{-})$.
Each of the tensors is sliced by the selected index, for instance,
$\boldsymbol{A}^{(1)}_{-}$ is sliced as $\boldsymbol{A}^{(1)}_{s_{1},-}$ for
fixed value $s_{1}\in\\{1,2,\ldots,L_{1}\\}$ of index $1$. Then, for fixed
value $s_{i}$ for each index $i$, the equation for each sliced tensor network
corresponding to the string $s_{1}s_{2}...s_{n}$ is

$\boldsymbol{\sigma}_{s_{1}s_{2}...s_{n}}=\sum_{-}(\prod_{j=1}^{m_{1}}\boldsymbol{A}^{(j)}_{s_{1},-}\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{s_{2},-}...\prod_{l=1}^{m_{n}}\boldsymbol{Z}^{(l)}_{s_{n},-})\boldsymbol{T}_{-},$

where each index $i$ spans by $s_{i}\in\\{1,2,\ldots,L_{i}\\}$. Hence, there
are $(N=\prod_{i=1}^{n}L_{i})$ product terms in the summation we need to
retrieve.

For an index $i$, we define a template polynomial $p_{i}(x)$ to encode tensors
sharing the corresponding index as

$p_{i}(x)=1+x+x^{1+m_{i}}+x^{1+m_{i}+m_{i}^{2}}+\ldots+x^{1+m_{i}+\dots+m_{i}^{L_{i}-2}}$
$=\sum_{j=1}^{L_{i}}x^{\frac{m_{i}^{j-1}-1}{m_{i}-1}}.$

For instance, encoding polynomial for index $1$ is $p_{1}(x)$. That is, we
encode each tensor sliced by the index $1$ by putting the sliced tensors on
each of the coefficients of encoding polynomial $p_{1}(x)$:

$\tilde{\boldsymbol{A}}^{(k)}_{-}=\sum_{j=1}^{L_{1}}\boldsymbol{A}^{(k)}_{j,-}x^{\frac{m_{1}^{j-1}-1}{m_{1}-1}}(k=1,2,..,m_{1})$
$=\boldsymbol{A}^{(k)}_{1,-}+\boldsymbol{A}^{(k)}_{2,-}x+\boldsymbol{A}^{(k)}_{3,-}x^{1+m_{i}}+\ldots+\boldsymbol{A}^{(k)}_{L_{1},-}x^{\frac{m_{1}^{L_{1}-1}-1}{m_{1}-1}},$

for $k\in\\{1,2,...,m_{1}\\}$.

Then, the product of the encoded tensors results in a polynomial as below:

$\prod_{k=1}^{m_{1}}\tilde{\boldsymbol{A}}^{(k)}=\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{1,-}+x^{m_{1}}\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{2,-}+x^{m_{1}^{2}+m_{1}}\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{3,-}$
$+x^{m^{3}_{1}+m^{2}_{1}+m_{1}}\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{4,-}+\ldots+x^{m_{1}^{L_{1}-1}+...+m_{1}}\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{L_{1},-}$
$+(\text{other order
terms})=poly(x^{m_{1}^{L_{1}-1}+...+m_{1}}=x^{\frac{m_{1}^{L_{1}}-m_{1}}{m_{1}-1}}).$

Each term of $\prod_{k=1}^{m_{1}}\boldsymbol{A}^{(k)}_{s_{1},-}$ for
$s_{1}\in\\{1,2,\ldots,L_{1}\\}$ is confined on different order term only by
itself. That is because order of individual term of the encoding polynomial
$p_{1}(x)$ is larger than $\times m_{1}$ of previous order.

For index $2$, we use encoding polynomial
$p_{2}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})$ where

$p_{2}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})=1+x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}+x^{(1+m_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}+x^{(1+m_{2}+m_{2}^{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}$
$+\ldots+x^{(1+m_{2}+\dots+m_{2}^{L_{2}-2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}=\sum_{j=1}^{L_{2}}x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\frac{m_{2}^{j-1}-1}{m_{2}-1}}.$

Same as for index $1$, we encode each tensor of index $2$ by putting the
sliced tensors on each of the coefficients of encoding polynomial
$p_{2}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})$:

$\tilde{\boldsymbol{B}}^{(k)}_{-}=\sum_{j=1}^{L_{2}}\boldsymbol{B}^{(k)}_{j,-}x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\frac{m_{2}^{j-1}-1}{m_{2}-1}}(k=1,2,..,m_{2})$
$=\boldsymbol{B}^{(k)}_{1,-}+\boldsymbol{B}^{(k)}_{2,-}x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}+\boldsymbol{B}^{(k)}_{3,-}x^{(1+m_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}$
$+\boldsymbol{B}^{(k)}_{4,-}x^{(1+m_{2}+m^{2}_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}+\ldots+\boldsymbol{B}^{(k)}_{L_{2},-}x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\frac{m_{2}^{L_{2}-1}-1}{m_{2}-1}}.$

Then the product of the encoded tensors for index $2$ results in a polynomial
as below:

$\prod_{k=1}^{m_{2}}\tilde{\boldsymbol{B}}^{(k)}=\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{1,-}+x^{m_{2}\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}\prod_{k=1}^{m_{1}}\boldsymbol{B}^{(k)}_{2,-}$
$+x^{(m_{2}^{2}+m_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}\prod_{k=1}^{m_{1}}\boldsymbol{B}^{(k)}_{3,-}+x^{(m^{3}_{2}+m^{2}_{2}+m_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}\prod_{k=1}^{m_{1}}\boldsymbol{B}^{(k)}_{4,-}$
$+\ldots+x^{(m_{2}^{L_{2}-1}+...+m_{2})\frac{m_{1}^{L_{1}}-1}{m_{1}-1}}\prod_{k=1}^{m_{1}}\boldsymbol{B}^{(k)}_{L_{2},-}+(\text{other
order terms})$
$=poly(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}\frac{m_{2}^{L_{2}}-m_{2}}{m_{2}-1}}).$

Just as the case of index $1$, each term of
$\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{s_{2},-}$ for
$s_{2}\in\\{1,2,\ldots,L_{2}\\}$ is confined on different order term only by
itself. Likewise, this is because order of individual term of the encoding
polynomial $p_{2}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})$ is larger than $\times
m_{2}$ of previous order.

Then, polynomial resulting from
$(\prod_{k=1}^{m_{1}}\tilde{\boldsymbol{A}}^{(k)}\prod_{k=1}^{m_{2}}\tilde{\boldsymbol{B}}^{(k)})$
confines every individual term of
$(\prod_{j=1}^{m_{1}}\boldsymbol{A}^{(j)}_{s_{1},-}\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{s_{2},-})$
for $s_{1}\in\\{1,2,\ldots,L_{1}\\},s_{2}\in\\{1,2,\ldots,L_{2}\\}$ each in
different order term. That is because order of individual terms of
$p_{2}(x^{\frac{m_{1}^{L_{1}}-1}{m_{1}-1}})$ is larger than
$[m_{2}\times(\text{previous order})+\frac{m_{1}^{L_{1}}-m_{1}}{m_{1}-1}]$.

To generalize, tensors sliced by index $i$ are encoded by encoding polynomial
of

$p_{i}(x^{\prod_{j=1}^{i-1}\frac{m_{j}^{L_{j}}-1}{m_{j}-1}}),$

by putting the sliced tensors on each of the coefficients of the encoding
polynomial. For instance, encoded tensors for last index $n$ are as follows:

$\tilde{\boldsymbol{Z}}^{(k)}_{-}=\sum_{j=1}^{L_{n}}\boldsymbol{Z}^{(k)}_{j,-}x^{\prod_{j=1}^{n-1}\frac{m_{j}^{L_{j}}-1}{m_{j}-1}}(k=1,2,..,m_{n})$

As a result, an encoded tensor network results as a polynomial at the end of
the contraction, which confines every desired outcome of sliced tensor
networks each in different order term as below:

$\tilde{\boldsymbol{\sigma}}=\sum_{-}(\prod_{j=1}^{m_{1}}\tilde{\boldsymbol{A}}^{(j)}_{s_{1},-}\prod_{k=1}^{m_{2}}\tilde{\boldsymbol{B}}^{(k)}_{s_{2},-}...\prod_{l=1}^{m_{n}}\tilde{\boldsymbol{Z}}^{(l)}_{s_{n},-})\boldsymbol{T}_{-}$
$=\sum_{-}[\prod_{j=1}^{m_{1}}\boldsymbol{A}^{(j)}_{1,-}\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{1,-}...\prod_{l=1}^{m_{n}}\boldsymbol{Z}^{(l)}_{1,-}+\ldots$
$+x^{(\prod_{j=1}^{n}\frac{m_{j}^{L_{j}}-1}{m_{j}-1})-1}\boldsymbol{A}^{(j)}_{L_{1},-}\prod_{k=1}^{m_{2}}\boldsymbol{B}^{(k)}_{L_{2},-}...\prod_{l=1}^{m_{n}}\boldsymbol{Z}^{(l)}_{L_{n},-}]\boldsymbol{T}_{-}$
$=\boldsymbol{\sigma}_{11...1}+\ldots+\boldsymbol{\sigma}_{L_{1}L_{2}...L_{n}}x^{(\prod_{j=1}^{n}\frac{m_{j}^{L_{j}}-1}{m_{j}-1})-1}$
$=poly(x^{(\prod_{j=1}^{n}\frac{m_{j}^{L_{j}}-1}{m_{j}-1})-1})$

Every coefficient of this resulting polynomial can be retrieved by polynomial
interpolation with $(\prod_{i=1}^{n}\frac{m_{i}^{L_{i}}-1}{m_{i}-1})$
different arbitrary constants. Hence, with any successful
$(\prod_{i=1}^{n}\frac{m_{i}^{L_{i}}-1}{m_{i}-1})$ encoded tensor networks
each using different arbitrary constant $x$, all the desired outcomes of
sliced tensor networks can be retrieved and the final outcome of original
tensor network can be acquired by their summation.

Therefore, for $f$ node failures, $f$-resilient number is

$f+\prod_{i=1}^{n}\frac{m_{i}^{L_{i}}-1}{m_{i}-1}.$

for hyperedge code.