



\section{Beyond Optimality}\label{sec:beyond}

The algorithm presented in Sec.~\ref{sec:optimal} yields optimal linear contraction orders for tree tensor networks. In contrast, optimal \textit{general} contraction ordering for tree tensor networks is still open. This is also true for join ordering, where finding the optimal bushy plan (the equivalent of general order) for tree queries is also still open. Even so, we can obtain near-optimal orders, as shown in~\cite{adaptive}, where the key idea is to use the linear solutions as a starting point for the cubic DP. The main motivation is shown in the same work: If we obtain the linearization that gives the best general contraction order, then the cubic DP will find the optimal solution. This is also the insight discovered later in~\cite{cubic_dp_in_tensor} in the context of tensor networks. The main advantage of this approach is that, although we cannot find the optimal \textit{general} contraction order, we are at least optimal in the regime of \textit{linear} contraction orders. Therefore, we can ensure the robustness of the generated contraction orders.

Additionally, we can extend the idea to arbitrary tensor networks by applying the heuristic presented in~\cite{ik}: Since the IKKBZ algorithm requires tree-shaped graphs, we can first find a suitable spanning tree according to a greedy strategy and finally use it as input to IKKBZ. The heuristic selects the minimum spanning tree since the join result has size $s_{i,j}|R_i||R_j|$, where $|R_i|$ and $|R_j|$ are the relations to be joined and $s_{i,j}$ is the join \textit{selectivity}. Hence, the goal is to minimize the output size. For the context of tensor networks, we propose the use of the \textit{maximum} spanning tree, since we can rewrite the contraction cost of two tensors $T_1 \in \mathbb{C}^{p_1\times \ldots \times p_a \times q_1 \times \ldots \times q_b}$ and $T_2 \in \mathbb{C}^{q_1 \ldots \times q_b \times r_1 \times \ldots \times r_c}$ as $ \frac{|T_1||T_2|}{\prod_{i=1}^{b} q_i}$, which allows us to interpret the product of the common legs as the inverse of the selectivity.
\section{Conclusion}\label{sec:conclusion}

Contraction orders play a fundamental role in efficient simulations of quantum many-body systems and quantum circuits. Since the problem of finding optimal contraction orders for arbitrary tensor networks is NP-complete, previous research has focused on finding heuristic algorithms that are experimentally shown to output good orders. In this work, we have taken a more conservative approach and shown that tree tensor networks accept optimal linear contraction orders. The key idea comes from a decades-old algorithm from database join ordering. We also propose two other heuristics for near-optimal contraction orders by taking advantage of advanced research in join ordering~\cite{adaptive} and adapting them to the setting of tensor networks. In this way, we bring together two separate lines of research, both of which have an open problem: Can we find optimal \textit{general} contraction orders for tree tensor networks?
\section{Results}\label{sec:results}

In this section, we highlight the advantage of the polynomial-time algorithm compared to a standard dynamic programming algorithm for tree tensor networks that outputs the optimal linear contraction order. In the following, we will refer to the IKKBZ adapted for tensor networks as IKS.

\begin{figure}
\centering\includegraphics[width=7.5cm,height=7.5cm,keepaspectratio]{photo/output-16.png}
\caption{Polynomial-time advantage over exponential-time DP for small-sized tensor networks.}
\label{fig:plot1}
\end{figure}

\begin{figure}
\centering\includegraphics[width=7.5cm,height=7.5cm,keepaspectratio]{photo/output-64.png}
\caption{Polynomial-time advantage over exponential-time DP for larger tree tensor networks. We set a timeout of 10 seconds for each algorithm.}
\label{fig:plot2}
\end{figure}

To this end, we generate tree tensor networks of size 5 up to 64, each with 100 instances, as presented in~\cite{adaptive}. We run DP and IKS for each size and report the average execution time. Since DP inherently takes exponential time, we set a timeout of 10 seconds. To show the benefits of the polynomial-time algorithm, we first plot the results for small tree-tensor networks in Fig.~\ref{fig:plot1} and then for larger tree-tensor networks in Fig.~\ref{fig:plot2}, where DP fails after a certain size due to the exponential nature of the problem.
\section{Introduction}\label{sec:introduction}
In this work, we consider tensors to be multi-dimensional arrays and define them as $T \in \mathbb{C}^{n_1 \times \ldots \times n_k}$. Tensor networks can be intuitively understood from the example of matrix multiplication: Consider two matrices $A\in\mathbb{C}^{p \times q}$ and $B\in\mathbb{C}^{q \times r}$, the product of which is defined as $C_{ij} = (AB)_{ij} = \sum_{k}A_{ik}B_{kj}$, where $C \in \mathbb{C}^{p \times r}$. Thus, there is a common index over which we perform the multiplication. In tensor notation, this represents the \textit{contraction} over their common leg (dimension) and can be represented in a graphical notation as in Fig.~\ref{fig:example}.

\input{figures/example.tex}

For two tensors $T_1 \in \mathbb{C}^{p_1 \times \ldots \times p_a \times q_1 \times \ldots \times q_b}$ and $T_2 \in \mathbb{C}^{q_1 \times \ldots \times q_b \times r_1 \times \ldots \times r_c}$, the result of the contraction is a tensor $T_{1,2}\in \mathbb{C}^{p_1 \times \ldots \times p_a \times r_1 \times \ldots \times r_c}$, i.e., the common legs, $q_i$, disappear, while the others, $p_i$ and $r_i$, respectively, remain.

To systematically represent contractions of a set of tensors, i.e., a tensor network, we define $G := (V, E, f)$ to be an undirected weighted graph, where $V$ is the set of tensors and $E$ the edge set, where $\{T_u, T_v\} \in E$ means that tensors $T_u$ and $T_v$ need to be contracted over the common leg of size $f({T_u, T_v})$. For simplicity, we only consider tensor networks whose output is a scalar. The proposed algorithms are not affected by this assumption, rather it saves us from additional notations.

To see why contraction order matters, consider the multiplication of three matrices $A\in\mathbb{C}^{20\times 30}$, $B\in\mathbb{C}^{30\times 10}$, and $C\in\mathbb{C}^{10 \times 50}$. There are two ways they can be multiplied: $(i)$ $(AB)C$ and $(ii)$ $A(BC)$. The cost of a matrix multiplication is defined as the number of scalar multiplications performed. For two arbitrary matrices $M_1\in\mathbb{C}^{p \times q}$ and $M_2\in\mathbb{C}^{q\times r}$ this is equal to $pqr$. Thus, variant $(i)$ incurs a cost of $20 \cdot 30 \cdot 10 + 20 \cdot 10 \cdot 50 = 16000$, while variant $(ii)$ incurs a cost of $30 \cdot 10 \cdot 50 + 20 \cdot 30 \cdot 50 = 45000$. We conclude that $(AB)C$ is the order to be preferred. In the same way, this can be generalized to tensor networks, where the cost of contracting tensors $T_1 \in \mathbb{C}^{p_1\times \ldots \times p_a \times q_1 \times \ldots \times q_b}$ and $T_2 \in \mathbb{C}^{q_1 \times \ldots \times q_b \times r_1 \times \ldots \times r_c}$ is equal to $p_1 \cdot \ldots \cdot p_a \cdot q_1 \cdot \ldots \cdot q_b \cdot r_1 \cdot \ldots \cdot r_c$. Hence, our goal is to find the \textit{optimal} contraction order, i.e., the one which minimizes the total cost.

We additionally differentiate between \textit{linear} contraction orders and \textit{general} ones, as established in~\cite{contraction_trees}. A linear contraction order contracts tensors one by one (think of a multiplication chain), while a generalized contraction order does not have any such restrictions. In the case of matrix multiplication of four matrices $A, B, C$, and $D$, an example of a linear contraction order is $(((AB)C)D)$, while an example of a general one is $(AB)(CD)$.
\newpage
The main drawback is that the problem of finding optimal contraction orders, either general or linear, for arbitrary tensor networks is NP-complete~\cite{Lam1997OnOA}. This means that we cannot ensure optimal orders for large tensors networks. For this reason, there exists a number of heuristics which promise good contraction orders.

We argue that this is a defeatist position and show that we can indeed find optimal contraction orders for a more restrictive class of tensor networks. We take a step back from employing heuristics and ask ourselves the following question: For which class of tensor networks and for what kind of contraction orders, i.e., linear or general, can we ensure optimality? A simple answer is to consider chain tensor networks, where a general contraction order can be found in polynomial time using the same cubic dynamic programming (DP) as for the matrix chain ordering problem. The answer we give in this work considers a broader class, namely \textit{tree} tensor networks, for which we show that we can guarantee optimal \textit{linear} contraction orders.

Beyond the optimality results, we highlight two join ordering techniques that can build on our work to guarantee near-optimal orderings for arbitrary tensor networks.

The rest of the paper is organized as follows: Section~\ref{sec:related} summarizes related work, both in tensor networks and database join ordering. In section~\ref{sec:optimal} we outline the algorithm for optimal linear contraction orders for tree tensor networks. In section~\ref{sec:beyond}, inspired from recent research on join ordering~\cite{adaptive}, we discuss algorithms for near-optimal contraction orders on arbitrary tensor networks. We summarize the results in section~\ref{sec:results}. Finally, we draw conclusions and discuss future work.
\section{Motivation}\label{sec:motivation}
Let us first motivate why we need tensor networks. They are approximations of high-order tensors, hence easier to handle. For this reason, they have recently found application in physics~\cite{Schollw_ck_2005}, quantum chemistry~\cite{Nakatani_2013}, and machine learning~\cite{ml_tensor}. In particular, for quantum simulations, tensor networks can effortlessly simulate a quantum circuit classically, i.e., without quantum computation. This is indeed promising since we are still in the NISQ era~\cite{Preskill_2018}, where a large number of qubits represents a challenge.

However, while the contraction operation itself is trivial, its execution time depends heavily on the order in which the contractions are executed. In the following sections, we introduce the concepts and formulate the problem.
\section{Optimal Contraction Orders}\label{sec:optimal}

As stated in the previous sections, we can indeed provide optimal \textit{linear} contraction orders for \textit{tree} tensor networks. This relies on the IKKBZ algorithm~\cite{ik,kbz} which itself builds upon the solution to a well-known job ordering problem~\cite{Monma1979SequencingWS}.

Let $T_{a,b} := T_aT_b$ denote the contraction of two tensors $T_a$ and $T_b$ and $c(T_aT_b)$ their contraction cost, as presented in Sec.~\ref{sec:introduction}. Further, we do not consider outer products, i.e., contracting tensors which do not share a common leg. In join ordering, these are called cross products, and the original IKKBZ algorithm cannot guarantee optimality in their presence.

The problem statement of optimal linear contraction ordering can then be stated as finding the permutation $\pi$ of tensors which minimizes the contraction cost, i.e.,

\begin{equation}
\argmin_{\pi} \Gamma(\pi) = \sum_{i=2}^{N} c(T_{\pi_1,\ldots,\pi_{i-1}}T_{\pi_i}).
\label{eqn:cost}
\end{equation}

Furthermore, as we deal with \textit{tree-shaped} tensor networks, it can be rewritten as 

\begin{equation}
\argmin_{\pi} \Gamma(\pi) = \sum_{i=2}^{N} \frac{|T_{\pi_1,\ldots,\pi_{i-1}}|| T_{\pi_i}|}{|e_{T_{\pi_j}, T_{\pi_i}}|},
\label{eqn:reformulated_cost}
\end{equation}

where $e_{T_{\pi_j}, T_{\pi_i}}$ is the \textit{unique} edge connecting $T_{\pi_i}$ to the already contracted tensor network $T_{\pi_1,\ldots,\pi_{j},\ldots,\pi_{i-1}}$. This uniqueness is implied by the fact that outer products are not allowed. We can ensure this constraint by introducing a partial order on the tensors, represented by a precedence graph:

\subsection{Precedence Graph}\label{sec:prec_graph}

One of the prerequisites for the optimality of the algorithm is the tree shape of the network. To achieve optimality, the algorithm roots the network in each tensor, solves the ordering problem for the obtained rooted tree, and picks the solution which leads to the minimum cost. The rooted tree is referred to as the precedence graph and describes the partial order of the contractions. An example for a precedence graph of a tree tensor network can be seen in Fig.~\ref{fig:precedence_graph}. In that example, the precedence graph specifies $T_4 \rightarrow T_3$ and $T_4 \rightarrow T_2$, i.e., $T_4$ should occur in the order \textit{before} $T_3$ and $T_2$, respectively.
\input{figures/precedence.tex}
\subsection{Rank Function}\label{subsec:rank}

The key idea behind the algorithm is that it allows global improvement based on local decisions. While this is not a novel idea in algorithmic design (see Kruskal's algorithm for the the minimum spanning tree~\cite{kruskal}), it is a powerful one. In our case, this is achieved through a $rank$ function. To illustrate the idea, let us consider the precedence graph depicted in Fig.~\ref{fig:example}, where $T_1\in\mathbb{C}^{p \times q}, T_2\in\mathbb{C}^{q \times r}, T_3\in\mathbb{C}^{r}$, and $T_4\in\mathbb{C}^{q}$.

\begin{figure}
\centering
\begin{tikzpicture}
\begin{scope}[auto, node distance=0.25cm,
   node_style/.style={align=center, circle,draw,inner sep=0pt, text centered,minimum size=0.5cm,
    font=\sffamily},
   edge_style/.style={draw=black}]


    \node[node_style,label={[xshift=0em, yshift=0.5em,font=\footnotesize,teal]}] (v1) at (1,5) {$T_1$};
    \node[node_style,label={[xshift=0em, yshift=0.5em,font=\footnotesize,teal]}] (v2) at (0.5,4.33) {$T_2$};
    \node[node_style,label={[xshift=0em, yshift=0.5em,font=\footnotesize,teal]}] (v3) at (1.5,4.33) {$T_3$};
    \node[node_style,label={[xshift=0em, yshift=-3.5em,font=\footnotesize,teal]}] (v4) at (0.5,3.5) {$T_4$};
    \draw[edge_style]  (v1) edge node{$p$} (v2);
    \draw[edge_style]  (v2) edge node{$q$} (v4);
    \draw[edge_style]  (v1) edge node{$r$} (v3);
    \end{scope}
\end{tikzpicture}
\caption[Rank Function]{Precedence graph of $T_1$, where $T_1\in\mathbb{C}^{p \times q}, T_2\in\mathbb{C}^{q \times r}, T_3\in\mathbb{C}^{r}$, and $T_4\in\mathbb{C}^{q}$}
\label{fig:example}
\end{figure}

Assume that we want to find out whether it is more advantageous to perform $(i)$ $T_1\begingroup\color{purple}T_2\endgroup \begingroup\color{cyan}T_4\endgroup T_3$ or $(ii)$ $T_1 
 \begingroup\color{cyan}T_4\endgroup \begingroup\color{purple}T_2\endgroup T_3$, i.e., which tensor between $\begingroup\color{purple}T_2\endgroup$ and $\begingroup\color{cyan}T_4\endgroup$ to contract first. Let us calculate the cost of each order:

\begin{enumerate}[$(i)$]
    \item $c(T_1T_{\begingroup\color{purple}2\endgroup})=pqr$, $c(T_{1,\begingroup\color{purple}2\endgroup}T_{\begingroup\color{cyan}4\endgroup})=qr$, and $c(T_{1,\begingroup\color{purple}2\endgroup,\begingroup\color{cyan}4\endgroup}T_3)=r$, yielding a total cost of $pqr + qr + r$,
    \item $c(T_1T_{\begingroup\color{cyan}4\endgroup})=pq$, $c(T_{1,{\begingroup\color{cyan}4\endgroup}}T_{\begingroup\color{purple}2\endgroup})=pr$, and $c(T_{1,{\begingroup\color{cyan}4\endgroup},{\begingroup\color{purple}2\endgroup}}T_3)=r$, yielding a total cost of $pq + pr + r$.
\end{enumerate}
When comparing both costs, i.e.,
\begin{flalign*}
    &\qquad& pqr + qr + r &\leq pq + pr + r &&\\
    \iff&& pqr + qr &\leq pq + pr &&\\
    \iff&& pq(r - 1) &\leq r(p - q), &&\\
   
\end{flalign*}
\label{eq:naive} we observe an interesting pattern: the original inequality reduces to a simpler one that considers only the tensors $\begingroup\color{purple}T_2\endgroup$ and $\begingroup\color{cyan}T_4\endgroup$ to be compared. To see this, assume $r > 1$ and $p > q$. Then the inequality is equivalent to $\frac{pq}{p - q} \leq \frac{r}{r - 1}$. Hence, a local decision is enough to decide whether the order $T_1\begingroup\color{purple}T_2\endgroup \begingroup\color{cyan}T_4\endgroup T_3$ is better than $T_1\begingroup\color{cyan}T_4\endgroup \begingroup\color{purple}T_2\endgroup T_3$.

Indeed, this pattern is not arbitrary: it is the core idea behind the Adjacent Sequence Interchange (ASI) property~\cite{Monma1979SequencingWS}, employed by the IKKBZ algorithm:

\begin{definition}[ASI Property]\label{def:asi}
Let $A$ and $B$ be two sequences and $U$ and $V$ two non-empty sequences. We say a cost function $C$ has the adjacent sequence interchange (ASI) property, iff there exists a function $rank$ such that the following holds
\[
C(AUVB) \leq C(AVUB) \iff rank(U) \leq rank(V),
\]
if AUVB and AVUB satisfy the precedence constraints imposed by a given precedence graph.
\end{definition}

Thus, a cost function $C$ with ASI property allows us to use \emph{local} comparisons based on the rank function, i.e., $rank(U) \leq rank(V)$, to infer comparisons between \emph{global} costs, i.e., $C(AUVB) \leq C(AVUB)$. Our algorithm will therefore make extensive use of rank comparisons to find optimal linear contraction orders.

Let $T_S$ be the tensor represented by the sequence $S$, where the leg pointing to the parent is not taken into account. Then, using the reformulation from Eq.~\ref{eqn:reformulated_cost}, we can prove that $\Gamma$, the total contraction cost, has the ASI property:

\begin{lemma}
$\Gamma$ has the ASI property.
\end{lemma}
\begin{proof}
{\small
\begin{align*}
\qquad&\Gamma(AUVB) \leq \Gamma(AVUB) &&\\    
\iff& \Gamma(AUV) + \frac{|T_{AUV}|}{|e_{AUV,B}|}\Gamma(B) \leq \Gamma(AVU) + \frac{|T_{AVU}|}{|e_{AVU,B}|}\Gamma(B) &&\\
\iff& \Gamma(AUV) \leq \Gamma(AVU) &&\\
\iff& \Gamma(AU) + \frac{|T_{AU}|}{|e_{AU,V}|}\Gamma(V) \leq \Gamma(AV) + \frac{|T_{AV}|}{|e_{AV,U}|}\Gamma(U) &&\\
\iff& \frac{|T_{A}|}{|e_{A,U}|}\Gamma(U) + \frac{|T_{AU}|}{|e_{AU,V}|}\Gamma(V) \leq \frac{|T_{A}|}{|e_{A,V}|}\Gamma(V) + \frac{|T_{AV}|}{|e_{AV,U}|}\Gamma(U) &&\\
\iff& \frac{\Gamma(U)}{|e_{A,U}|} + \frac{|T_{U}|}{|e_{A\begingroup\color{red}U\endgroup,V}||e_{A,U}|}\Gamma(V) \leq \frac{\Gamma(V)}{|e_{A,V}|} + \frac{|T_{V}|}{|e_{A\begingroup\color{red}V\endgroup,U}||e_{A,V}|}\Gamma(U) &&\\
\iff& \frac{\Gamma(U)}{|e_{A,U}|} + \frac{|T_{U}|}{|e_{A,V}||e_{A,U}|}\Gamma(V) \leq \frac{\Gamma(V)}{|e_{A,V}|} + \frac{|T_{V}|}{|e_{A,U}||e_{A,V}|}\Gamma(U) &&\\
\iff& \frac{\Gamma(U)}{|e_{A,U}|}(1 - \frac{|T_V|}{|e_{A,V}|}) \leq \frac{\Gamma(V)}{|e_{A,V}|}(1 - \frac{|T_{V}|}{|e_{A,U}|}) &&\\
\iff& \Gamma(U)(|e_{A,V}| - |T_V|) \leq \Gamma(V)(|e_{A,U}| - |T_{V}|) &&\\
\end{align*}
}
We can in principle define $rank(S) := \frac{\Gamma(S)}{|e_{S}| - |T_S|}$, but this is not allowed since in case the denominator is negative, the sign of the inequality is reversed. Therefore, we work directly with the fraction representing the $rank$ and let $d_S := \Gamma(S)$ and $n_S := |e_{S}| - |T_S|$ be the denominator and numerator, respectively. We will thus evaluate  $rank(U) \leq rank(V)$ as $d_U \cdot n_V \leq d_V \cdot n_U$ so that we can directly work with the original inequality.
\end{proof}

\subsection{Algorithm}

As mentioned in Sec.~\ref{sec:related}, we adapt the IKKBZ algorithm to the setting of tree tensor networks. This has already been achieved by defining the corresponding $rank$ function. Apart from that, the algorithm remains the same as formulated in~\cite{ik,kbz}. In the following, the algorithm is outlined.

The algorithm treats each tensor $T_i \in V$ separately, by executing two consecutive phases:

\textbf{First Phase:} Build the precedence graph of $T_i$, as described in Sec.~\ref{subsec:rank}.\par

\textbf{Second Phase:} Compute the optimal permutation, which reflects the order in which the tensors should be contracted. The process of transforming the precedence graph into a permutation is called linearization and can be regarded as a recursive method: once the subtrees of the children are linearized, we can obtain the linearization of the whole tree by interleaving (merging) the tensors according to their ranks.

To this end, the notion of contradictory or conflicting sequences comes into play: whenever the precedence graph requires $A \rightarrow B$, but $rank(A) \geq rank(B)$, we combine $A$ and $B$ into a new single node. This node represents a \textit{compound} tensor which includes all tensors in $A$ and $B$. This operation, called normalization, is applied after interleaving the linearizations and preserves the invariant that ranks of a linearization are in \textit{increasing} order (cf. Def.~\ref{def:asi}). Thus, the only contradictory sequences we might obtain are between the root of a subtree and the first (compound) tensors in the result of the merge step. In that case, several normalization steps are applied, until no contradictory sequence is due.

To see why normalization steps are required for optimality, consider the following Lemma:
\begin{lemma}\label{lemma:normalize}
Let $\{A,B\}$ be a compound tensor. If $A \rightarrow B$ and $rank(A) \geq rank(B)$, then we find an optimal sequence in which $B$ directly follows $A$.
\end{lemma}

\begin{proof}
Assume, by contradiction, that $B$ does not directly follow $A$. Since $A \rightarrow B$, any optimal sequence must be of the form $UAVBW$, with $V \neq \epsilon$.

\textit{Case} $rank(A) \geq rank(V)$: We can exchange $V$ and $A$ without increasing the costs.

\textit{Case} $rank(A) \leq rank(V)$: Since $rank(B) \leq rank(A)$, due to transitivity we have $rank(B) \leq rank(V)$. Thus, we can exchange $B$ and $V$ without increasing the costs. \qedhere
\end{proof}

At the end, the reverse operation, called denormalization, is performed on the linearization of the precedence graph, where each compound tensor is replaced by the sequence of tensors it contains. The obtained sequence of tensors represents exactly the optimal permutation for the given precedence graph.\par

Finally, after applying both phases for each tensor, the linearization with the lowest cost is the optimal solution.

As regards time complexity, the algorithm builds the solution in polynomial time. The first phase, in which the precedence graph is constructed, can be completed via Depth First Search (DFS) in linear time, starting from the root and directing the edges towards the unvisited nodes of the tree. In the second phase, while normalization is a cheap operation, requiring only linear time, linearization is the most expensive one. Interleaving $m$ linearizations of size $s$ each requires $\mathcal{O}(ms\log_2(m))$ operations. This can be done with a min-heap data structure of size $m$. The worst-case time complexity of each iteration is $\mathcal{O}(n\log_2(n))$, so $\mathcal{O}(n^2\log_2(n))$ in total, since the algorithm repeats the two phases for each of the $n$ tensors.

Later in~\cite{kbz}, the algorithm was also improved, its time complexity being reduced to $\mathcal{O}(n^2)$, by observing that linearizations of neighboring nodes are similar.

\subsection{Example}

\input{figures/IKKBZ.tex}

A sample algorithm execution can be visualized in Fig.~\ref{fig:iks_exec}. We are considering only the precedence graph of $T_4$, but the algorithm would have to consider all precedence graphs. In the first phase, the algorithm builds the precedence graph of $T_4$, which is equivalent to rooting the network in $T_4$. The subtree of $T_2$ is not yet linearized, thus the algorithm merges $T_1$ and $T_5$, by observing the rank monotonicity, i.e., $\cg{violet}{2}\:\cdot\:\cg{brown}{0} \leq \cg{violet}{1}\:\cdot\:\cg{brown}{2}$. After the merge step, a normalization step is required (step \fbox{\scriptsize 3}), since $T_2 \rightarrow T_1$, but $rank(T_2) \geq rank(T_1)$. The algorithm is then able to linearize $T_4$ and denormalize it (the compound tensor $T_2T_1$ is decomposed into individual tensors, as in step \fbox{\scriptsize 5}). Finally, we can build the linear contraction order from the denormalized linearization.
\section{Related Work}\label{sec:related}

As mentioned in the introduction, the problem of finding optimal contraction orders for arbitrary tensor networks is NP-complete. Therefore, one must settle for a mixture of exponential algorithms for small problems, e.g., $n \leq 20$, and otherwise hope for good contraction orders. For this reason, previous research has focused on the latter part, trying to find better heuristics. However, there has also been improvements for small-sized instances, for which novel dynamic programming techniques have been developed~\cite{Pfeifer_2014}.

In contrast, database join ordering has already benefited from these algorithms decades ago. Selinger presented the first dynamic programming algorithm for join ordering as early as 1970~\cite{selinger}. Later, the IKKBZ algorithm~\cite{ik,kbz}, which we adjust for tensor networks, adapted a well-known problem from job ordering~\cite{Monma1979SequencingWS} to provide optimal \textit{linear} plans for \textit{tree} queries. Similarly, we extend it to tensor networks to provide the same optimality guarantees for \textit{tree} tensor networks and \textit{linear} contraction orders. The challenge here, as it was decades ago with the IKKBZ algorithm, is to adapt the mathematical foundations to our setting.

In the field of tensor networks, we have seen DP algorithms that also emulate graph structure~\cite{Pfeifer_2014}; in join ordering there is~\cite{new_dp}. Several heuristics have been proposed~\cite{Gray_2021}, as was the case in join ordering~\cite{goo}. Interestingly, both research directions have recently used learned models to predict the orders~\cite{Marcus_2018,rl_tensor}. The work closest to ours is~\cite{ttn_poly}, which \textit{relaxes} the problem by using a different cost function: minimizing the most expensive contraction cost. In this way, one can actually find a polynomial-time algorithm for tree tensor networks. In contrast, our algorithm works directly with the original cost function.

We are not the first work to recognize the similarities between join ordering and contraction ordering. Dudek et al.~\cite{tensor_join_connection} have already highlighted this and used tree decompositions to solve the original problem, still relying on heuristic algorithms.

In the next sections we provide the modifications needed to enable the IKKBZ algorithm to work on tree tensor networks.
