ReLU Neural Networks of Polynomial Size for Exact Maximum Flow Computation \thanks{We thank Martin Skutella and Arturo Merino for many inspiring and fruitful discussions. We thank L{\'a
\begin{minipage
This paper studies the expressive power of artificial neural networks (NNs) with rectified linear units. To study them as a model of real-valued computation, we introduce the concept of Max-Affine Arithmetic Programs and show equivalence between them and NNs concerning natural complexity measures. We then use this result to show that two fundamental combinatorial optimization problems can be solved with polynomial-size NNs, which is equivalent to the existence of very special strongly polynomial time algorithms. First, we show that for any undirected graph with $n$ nodes, there is an NN of size $\mathcal{O(n^3)$ that takes the edge weights as input and computes the value of a minimum spanning tree of the graph. Second, we show that for any directed graph with $n$ nodes and $m$ arcs, there is an NN of size $\mathcal{O(m^2n^2)$ that takes the arc capacities as input and computes a maximum flow. These results imply in particular that the solutions of the corresponding parametric optimization problems where all edge weights or arc capacities are free parameters can be encoded in polynomial space and evaluated in polynomial time, and that such an encoding is provided by an NN.
neuralmaxflow.tex
\documentclass[11pt, a4paper]{article}

\title{ReLU Neural Networks of Polynomial Size for Exact Maximum Flow Computation
\thanks{We thank Martin Skutella and Arturo Merino for many inspiring and fruitful discussions. We thank L{\'a}szl{\'o} V{\'e}gh, Max Klimm, and Jennifer Manke for many valuable comments.}}


\author{
	\begin{minipage}{0.4\textwidth}
		\centering
		Christoph Hertrich
		\thanks{Funded by DFG-GRK 2434 Facets of Complexity.}\\
		Technische Universit\"at Berlin\\ \texttt{hertrich@math.tu-berlin.de}
	\end{minipage}\qquad
	\begin{minipage}{0.4\textwidth}
		\centering
		Leon Sering
		\thanks{Funded by DFG Excellence Cluster MATH+ (EXC-2046/1, project ID: 390685689). Most of the results in this paper have been developed while Leon Sering was employed at TU Berlin.}\\
		ETH Z\"urich\\ \texttt{leon@sering.eu}
	\end{minipage}
}
\date{}

\usepackage{hyperref}
\usepackage{bm, braket, amsmath, amssymb, amsthm, thm-restate, mathtools, enumerate}
\usepackage{cleveref}
\allowdisplaybreaks

\newcommand{\mypara}[1]{\subsubsection*{#1}}

\bibliographystyle{plainurl}

\usepackage[top=3.4cm, bottom=4.3cm]{geometry}
\usepackage[utf8]{inputenc}
\usepackage[linesnumbered, ruled, lined, boxed, noend]{algorithm2e}

\newcommand{\R}{\mathbb{R}}
\newcommand{\C}{\mathbb{C}}
\newcommand{\N}{\mathbb{N}}
\newcommand{\No}{\N_0}
\newcommand{\Q}{\mathbb{Q}}
\newcommand{\Z}{\mathbb{Z}}
\newcommand*\diff{\mathop{}\!\mathrm{d}}
%\renewcommand{\P}{\mathcal{P}}
\newcommand{\E}{\mathbb{E}}

\DeclareMathOperator{\dist}{dist}
\DeclareMathOperator{\MST}{MST}

\newcommand{\bigO}{\mathcal{O}}
\newcommand{\floor}[1]{\left \lfloor #1 \right \rfloor }
\newcommand{\abs}[1]{\left\lvert #1 \right\rvert}
\newcommand{\NP}{\ensuremath{\mathsf{NP}}\xspace}
\renewcommand{\P}{\ensuremath{\mathsf{P}}\xspace}
\newcommand{\OPT}{\ensuremath{\mathsf{OPT}}\xspace}

\renewcommand{\epsilon}{\varepsilon}

\newcommand{\cupdot}{\mathbin{\mathaccent\cdot\cup}}

\renewcommand{\l}{\ell}


\newtheorem{theorem}{Theorem}
\newtheorem{proposition}[theorem]{Proposition}
\newtheorem{corollary}[theorem]{Corollary}
\newtheorem{claim}[theorem]{Claim}

\Crefname{algocf}{Algorithm}{Algorithms}

\newenvironment{subproof}[1][\proofname]{%
	\renewcommand{\qedsymbol}{$\blacksquare$}%
	\begin{proof}[#1]%
	}{%
	\end{proof}%
}

\begin{document}

\maketitle

\begin{abstract}
	This paper studies the expressive power of artificial neural networks (NNs) with rectified linear units. To study them as a model of real-valued computation, we introduce the concept of \emph{Max-Affine Arithmetic Programs} and show equivalence between them and NNs concerning natural complexity measures. We then use this result to show that two fundamental combinatorial optimization problems can be solved with polynomial-size NNs, which is equivalent to the existence of very special strongly polynomial time algorithms. First, we show that for any undirected graph with $n$ nodes, there is an NN of size $\mathcal{O}(n^3)$ that takes the edge weights as input and computes the value of a minimum spanning tree of the graph. Second, we show that for any directed graph with $n$ nodes and $m$ arcs, there is an NN of size $\mathcal{O}(m^2n^2)$ that takes the arc capacities as input and computes a maximum flow. These results imply in particular that the solutions of the corresponding parametric optimization problems where all edge weights or arc capacities are free parameters can be encoded in polynomial space and evaluated in polynomial time, and that such an encoding is provided by an NN.
\end{abstract}

\section{Introduction}
Artificial neural networks (NNs) achieved breakthrough results in various application domains like computer vision, natural language processing, autonomous driving, and many more~\cite{LeCunBengioHinton:DeepLearning}.
Also in the field of combinatorial optimization (CO), promising approaches to utilize NNs for problem solving or improving classical solution methods have been introduced~\cite{BengioLodiProuvost:MLforCO}.
However, the theoretical understanding of NNs still lags far behind these empirical successes.

When talking about NNs, we always refer to \emph{feedforward NNs with rectified linear unit (ReLU) activations}, one of the most popular models in practice~\cite{glorot2011deep}. These NNs are directed, acyclic, computational graphs where each node computes an affine transformation of the outputs of its predecessors and applies the ReLU activation function $x\mapsto \max\{0,x\}$ on top. The full NN then computes a function mapping real-valued inputs to real-valued outputs. An important theoretical question about these NNs is concerned with their expressivity: which functions can be (exactly) represented by an NN of a certain size?

In order to tackle this question, we view NNs as a model of computation operating on real numbers. This is similar to Boolean circuits, or, even closer, arithmetic circuits, which are well-studied objects in complexity theory. We then investigate the computational complexity of different CO problems within this model.
Besides understanding the expressivity of NNs, we believe that this point of view is of high interest for the CO community, because it provides new insights on structural properties of classical CO problems and, eventually, on the interplay between NNs and traditional algorithms.

\mypara{Computational Power of Neural Networks}
Arora et al.\ \cite{Arora:DNNwithReLU} proved that the set of functions computable by ReLU~NNs (of any size) is precisely the set of continuous, piecewise linear (CPWL) functions. Even though their construction shows that logarithmic depth (in the input dimension) is always sufficient, these networks might have huge width such that no bounds on the total size can be obtained.
In particular, it remains an open question to identify classes of CPWL functions that are representable with polynomial-size NNs. In this paper, we study this question for functions corresponding to the solution of CO problems.

If there are polynomial-size NNs to solve a certain problem, then there exists a strongly polynomial time algorithm for that problem, simply by executing the NN. However, the converse might be false. This is due to the fact that ReLU NNs only allow a very limited set of possible operations, namely affine combinations and maxima computations. In particular, any function computed by such NNs is continuous, making it impossible to realize instructions like a simple \textbf{if}-branching based on a comparison of real numbers. In fact, there are models of computation for which the use of branchings is exponentially powerful~\cite{jerrum1982some}.

For some CO problems, classical algorithms do not involve conditional branchings and, thus, can easily be implemented as an NN. This is, for example, true for many dynamic programs. In these cases, the existence of efficient NNs follows immediately. We refer to Hertrich and Skutella~\cite{knapsackPaper} for some examples of this kind. In particular, polynomial-size NNs to compute the length of a shortest path in a network from given arc lengths are possible.

For other problems, like the Minimum Spanning Tree Problem or the Maximum Flow Problem, all classical algorithms seem to use conditional branchings.
For example, many maximum flow algorithms use them to decide whether an arc is part of the \emph{residual network}. Thus, these algorithms cannot be implemented with NNs and, without the results of this paper, it remains unclear whether polynomial-size NNs to solve these problems do exist at all.

More generally, the grand vision to which we contribute in this paper is to understand for different CO problems whether they can be solved by a limited class of strongly polynomial time algorithms that corresponds to polynomial-size NNs.

\mypara{Exact Neural Networks for CO Problems}
Many functions related to CO~problems are CPWL and can thus be represented with an NN (of any size), using the result by Arora et al.~\cite{Arora:DNNwithReLU}. In fact, for a feasible set $\mathcal{X}\subseteq\set{0,1}^n$, consider the generic CO~problem $\min_{x\in\mathcal{X}} c^T x$. This formulation covers a broad range of CO~problems, among them efficiently solvable problems like Shortest Path and Matching Problems, the Minimum Spanning Tree Problem, the Maximum Flow Problem (via duality to the Minimum Cut Problem), and also NP-hard problems like the Knapsack or the Traveling Salesperson Problem. When considering the cost vector $c$ to be a variable, the function $c\mapsto\min_{x\in\mathcal{X}} c^T x$ mapping~$c$ to the objective value, is a minimum of finitely many linear functions and thus CPWL.
For NP-hard problems it is unlikely that polynomial-size NNs can compute this function. However, even for problems in P, the existence of polynomial-size NNs is unclear as argued above.

Note that, due to their continuous nature, ReLU NNs (without threshold gates or similar) that output the discrete solution vector $x\in\set{0,1}^n$, and not only the objective value, cannot exist. Therefore, we are mostly interested in computing only the objective value. However, for the Maximum Flow Problem the situation is slightly different since it is not the primal Maximum Flow Problem but the dual Minimum Cut Problem that is represented in the generic form above. 
In fact, using an appropriate tie-breaking rule between equally good flows, the vector denoting a maximum flow is a CPWL function of the arc capacities.
Hence, a ReLU NN computing this maximum flow function of type~$\R^{m}\to\R^{m}$ ($m$~being the number of arcs in the flow-network) is possible, in addition to only computing the objective value.
In this paper, we show that this is even possible with polynomial size.

\mypara{Contributions and Outline}

In \Cref{Sec:RelWork}, we relate NNs to other models in complexity theory and survey further related work. In \Cref{Sec:Prelim}, we formally introduce ReLU~NNs. Afterwards, we provide the following specific contributions.

\begin{itemize}
	\item In order to make it possible to think about NNs in an algorithmic way, we introduce the pseudo-code language \emph{Max-Affine Arithmetic Programs} (MAAPs) in \Cref{Sec:MAAP}. We show that MAAPs and NNs are basically equivalent (up to constant factors) concerning three basic complexity measures corresponding to depth, width, and overall size of NNs. Hence, MAAPs serve as a convenient tool for constructing NNs with bounded size and could be useful for further research about NN expressivity far beyond the scope of this paper.
	\item Based on a result from subtraction-free circuit complexity~\cite{fomin2016subtraction}, we prove that polynomial-size NNs can exactly compute the value of a minimum spanning tree in \Cref{Sec:CO}. For a graph with~$n$ nodes, the size of the NN is~$\mathcal{O}(n^3)$.
	\item We prove that polynomial-size NNs can exactly solve the Maximum Flow Problem in \Cref{Sec:MaxFlow}. Since classical maximum flow algorithms do not admit straightforward NN implementations, we develop a new maximum flow algorithm in the form of a MAAP. In particular, our algorithm only consists of simple arithmetic computations and works without conditional branchings. We conclude that, for a flow network with $n$ nodes and $m$ arcs, there exists a ReLU NN of size $\mathcal{O}(n^2 m^2)$ that maps any real-valued arc capacities to a maximum flow.
\end{itemize}

We would like to mention that the latter result also has an interesting interpretation from the perspective of parametric algorithms. There exists a variety of literature concerning the question how one could represent solutions to the Maximum Flow Problem if the input depends on one or several unknown parameters; see, e.g., the works by Gallo~et~al.~\cite{gallo1989fast} and McCormick~\cite{mccormick1999fast}. Our results imply that there exists such a representation of polynomial size for the most general form of parametric maximum flow problems, namely the one where \emph{all} arc capacities are independent free parameters. This representation is given in the form of an NN and can be evaluated in polynomial time.

\section{Related Models and Literature}\label{Sec:RelWork}

\mypara{Boolean Circuits}
Even though NNs are naturally a model of real computation, it is worth to have a look at their computational power with respect to Boolean inputs. It is easy to see that ReLU~NNs can directly simulate \mbox{AND-,} \mbox{OR-,} and NOT-gates, and thus also any Boolean circuit; see~\cite{mukherjee2017lower}. Hence, in Boolean arithmetics, any problem in P can be solved with polynomial-size NNs. However, to grasp the full nature of NNs, one needs to consider real instead of Boolean arithmetics. Solving a problem for all real-valued inputs is a much stronger requirement and, consequently, the class of functions representable with polynomial-size NNs is much less understood than in Boolean arithmetics.

\mypara{Neural Networks and Arithmetic Circuits}
In a broad sense, \emph{arithmetic circuits} are directed acyclic graphs where each node computes some expression from the outputs of all its predecessors. Then, the full graph defines a function mapping the inputs, located at the sources of the graph, to the outputs, located at the sinks of the graph. In ordinary arithmetic circuits, the inner nodes can only compute either the sum or the product of the outputs of their predecessors. There has been a lot of research about the complexity of such arithmetic circuits~\cite{shpilka2010arithmetic}.
Closer to ReLU NNs, there is a special kind of arithmetic circuits called \emph{tropical circuits}~\cite{jukna2015lower}. In contrast to ordinary arithmetic circuits, they only contain maximum (or minimum) gates instead of sum gates and sum gates instead of product gates. Thus, they are arithmetic circuits in the max-plus algebra.

A tropical circuit can be simulated by an NN of roughly the same size since NNs can compute maxima and sums. Thus, NNs are at least as powerful as tropical circuits.
However, NNs are strictly more powerful. In particular, lower bounds on the size of tropical circuits do not apply to NNs. A particular example is the computation of the value of a minimum spanning tree. By Jukna and Seiwert~\cite{jukna2019greedy}, no polynomial-size tropical circuit can do this. However, as we show in \Cref{Sec:CO} as a consequence of a result by Fomin et al.~\cite{fomin2016subtraction}, an NN of cubic size (in the number of nodes of the input graph) is sufficient for this task.

The reason for this exponential gap is that, by using negative weights, NNs can realize subtractions (that is, tropical division), which is not possible with tropical circuits; compare the discussion by Jukna and Seiwert~\cite{jukna2019greedy}. However, this is not the only feature that makes NNs more powerful than tropical circuits. In addition, NNs can realize scalar multiplication (tropical exponentiation) with arbitrary real numbers via their weights, which is impossible with tropical circuits. It is unclear to what extend this feature increases the computational power of NNs compared to tropical circuits.

For these reasons, lower bounds from arithmetic circuit complexity do unfortunately not transfer to NNs. Hence, we identify it as a major challenge to prove meaningful lower bounds of any kind for the computational model of NNs.


\mypara{Neural Networks and Parallel Computation}
Similar to Boolean circuits, NNs can be seen as a model of parallel computation since the execution of all neurons of one layer can be done in parallel. Without going into detail here, the depth of an NN is related to the running time of a parallel algorithm, its width is related to the required number of processing units, and its size to the total amount of work conducted by the algorithm. Against this background, a natural goal is to design NNs as shallow as possible in order to make maximal use of parallelization. However, several results in the area of NN expressivity state that decreasing the depth is often only possible at the cost of an \emph{exponential} increase in width; see~\cite{Arora:DNNwithReLU,eldan2016power,liang2017deep,safran2017depth,Telgarsky15,telgarsky2016benefits,yarotsky2017error}.

Interestingly, a related observation can be made for the Maximum Flow Problem using complexity theory. 
By the result of Arora et al.\ \cite{Arora:DNNwithReLU} mentioned above, any CPWL function can be represented with logarithmic depth in the input dimension. In particular, this is also true for the function mapping arc capacities to a maximum flow. Hence, NNs with logarithmic depth can solve the Maximum Flow Problem and it arises the question whether such shallow NNs are also possible while maintaining polynomial total size.

The answer is most likely ``no'' since the Maximum Flow Problem is \mbox{\emph{P-complete}}~\cite{goldschlager1982}. P-complete problems are those problems in P that are \emph{inherently sequential}, meaning that there cannot exist a parallel algorithm with polylogarithmic running time using a polynomial number of processors unless the complexity classes P and NC coincide, which is conjectured to be not the case~\cite{greenlaw1995limits}. NNs with polylogarithmic depth and polynomial total size that solve the Maximum Flow Problem, however, would translate to such an algorithm (under mild additional conditions, such as, that the weights of the NN can be computed in polynomial time).
Therefore, we conclude that it is unlikely to obtain NNs for the Maximum Flow Problem that make significant use of parallelization.
In other words, NNs with polylogarithmic depth and polynomial width probably do not exist.

\mypara{Further Related Work}

NNs have been studied from a circuit complexity point of view before~\cite{beiu1996circuit,parberry1994circuit,shawe1992classes}. However, these works focus on Boolean circuit complexity of NNs with sigmoid or threshold activation functions. We are not aware of previous work investigating the computational power of ReLU~NNs as arithmetic circuits operating on the real numbers.

Using NNs to solve CO problems started with so-called \emph{Hopfield networks}~\cite{hopfield1985neural} and related architectures in the 1980s and has been extended to general nonlinear programming problems later on~\cite{kennedy1988neural}. Smith~\cite{Smith:NNforCOreview} surveys these early approaches. Also, specific NNs to solve the Maximum Flow Problem have been developed before~\cite{ali1991neural,effati2008neural,nazemi2012}. However, the NNs used in these works are conceptually very different from modern feedforward NNs that are considered in this paper.

In recent years interactions between NNs and CO have regained a lot of attention in the literature~\cite{BengioLodiProuvost:MLforCO}. For example, NNs have successfully been used to decide on which variables to branch when solving mixed-integer linear programs~\cite{lodi2017learning}. Also, a large variety of NN approaches for specific CO problems has been proposed~\cite{Bello:NeuralCOwithRL,Emami:learningPermutations,Khalil:LearningCOoverGraphs,kool2019attention,nowak:quadrAssignment,Vinyals:PointerNetworks}. However, these approaches have in common that they are of heuristic nature and no guarantees on the solution quality can be given.

Concerning the general expressivity of NNs, so-called \emph{universal approximation theorems}~\mbox{\cite{anthony2009neural,cybenko1989approximation,hornik1991approximation}} state that one hidden layer is already sufficient to approximate any continuous function on a compact domain arbitrarily well. Unfortunately, no insights concerning exact representability can be obtained from that. Recent progress concerning the classes of CPWL functions representable by ReLU neural networks with a certain architecture was made by Hertrich et al.~\cite{hertrich2021towards}. Various trade-offs between depth and width of NNs~\cite{Arora:DNNwithReLU,eldan2016power,hanin2019universal,hanin2017approximating,liang2017deep,nguyen2018neural,raghu2017expressive,safran2017depth,Telgarsky15,telgarsky2016benefits,yarotsky2017error} and approaches to count and bound the number of linear regions of a ReLU NN~\cite{hanin2019complexity,montufar2014regions,pascanu2014number,raghu2017expressive,serra2018bounding} have been found.

For an introduction to classical maximum flow algorithms, we refer to standard textbooks~\mbox{\cite{AhujaMagnantiOrlin:NetworkFlows,kortevygen,williamson_2019}}. The former two of these books also cover minimum spanning tree algorithms. The asymptotically fastest known maximum flow algorithm runs in $\mathcal{O}(nm)$ time for networks with $n$ nodes and~$m$ arcs~\cite{orlin2013}. However, as outlined above, classical algorithms are not applicable for direct implementation on an NN because they require conditional branchings. Nevertheless, our approach uses ideas of the famous algorithms by Edmonds-Karp~\cite{edmonds1972theoretical} and Dinic~\cite{dinic1970algorithm}.
	
	\section{Neural Networks with Rectified Linear Units}
	\label{Sec:Prelim}
	In this section, we formally define the primary object of study in this paper, using notations similar to \cite[Chapter 20]{Shalev2014:UnderstandingML}. A \emph{feedforward neural network with rectified linear units}, abbreviated by ReLU NN or simply NN, is a directed acyclic graph~$(V,E)$, for which each arc $e \in E$ is equipped with a weight $w_e\in\R$ and each
	node~$v \in V\setminus V_0$ is equipped with a bias $b_v\in\R$. Here, $V_0$ denotes the set of all nodes with no incoming arcs. The nodes~$V$ of an NN are
	called \emph{neurons} and the \emph{depth}~$k$ is given by the length of a longest path. The
	neurons are partitioned into \emph{layers} $V=V_0\cupdot V_1 \cupdot \cdots \cupdot V_k$ in such a way that the layer
	index strictly increases along each arc.\footnote{In other literature arcs are only allowed between successive layers.
		Clearly, this can always be achieved by introducing additional neurons. For our purposes, however, we want to avoid
		this restriction.} In addition, we assume that~$V_k$ is exactly the set of neurons with out-degree zero. The neurons in
	$V_0$ and $V_k$ are called \emph{input neurons} and \emph{output neurons}, respectively.
	All other neurons in $V\setminus(V_0\cup V_k)$ are so-called \emph{hidden neurons}. Let~$n_\l\coloneqq \abs{V_\l}$
	be the number of neurons in layer $\l$. The \emph{width} and \emph{size} of the NN are given by
	$\max\set{n_1,\dots,n_{k-1}}$ and~$\sum_{\l=1}^{k-1} n_\l$, respectively.
	
	In our paper it is crucial to distinguish fixed parameters of NNs, like architectural details and weights, from input, activation, and output values of neurons. We denote the latter by bold symbols in order to make the difference visible.
	
	The \emph{forward pass} of an NN is given by a function $\R^{n_0}\to \R^{n_k}$ that is obtained as follows. For an input vector
	$\bm x\in\R^{n_0}$ we inductively compute an \emph{activation} $\bm a(v)$ for every~\mbox{$v\in V\setminus V_0$} and an \emph{output} $\bm o(v)$
	for every \mbox{$v\in V\setminus V_k$}. First, the output values $\bm o(v)$ of the input neurons are set to the corresponding entries of the input vector $\bm x$. Second, the activation of a neuron~$v\in
	V\setminus V_0$ is given by the weighted sum of the outputs of all of its predecessors plus the bias $b_v$. Formally, we set $\bm a(v)\coloneqq b_v +
	\sum_{u\colon uv\in E} w_{uv} \bm o(u)$. Third, we apply the so-called \emph{activation function} $\sigma$ to obtain the output of all hidden neurons, i.e., $\bm o(v)\coloneqq\sigma(\bm a(v))$. In this work we only consider the \emph{ReLU function} $\sigma(z)=\max\set{0,z}$ as activation function. Finally, the activation values $\bm a(v)$ of the output neurons $v\in V_k$ provide the \emph{output vector} $\bm y\in \R^{n_k}$. Note that the activation function is not applied to these last activation values.
	An example of a very simple NN is given in \Cref{Fig:Min2Num}.
	%
	\begin{figure}[t]
		\centering
		\includegraphics[width=.25\linewidth]{minNN}
		%	\begin{tikzpicture}[scale=\tikzscalefactor, every node/.style={transform shape}]\tikzsize
			%		\node[smallneuron] (x1) at (0,5ex) {$\bm x_1$};
			%		\node[smallneuron] (x2) at (0,0) {$\bm x_2$};
			%		\node[smallneuron] (n1) at (11ex,5ex) {\relu};
			%		\node[smallneuron] (y) at (22ex,0ex) {$\bm y$};
			%		\scriptsize
			%		\draw[connection] (x1) -- node[above, pos=0.4]{-1} (n1);
			%		\draw[connection] (x2) -- node[above, pos=0.4]{1} (n1);
			%		\draw[connection] (n1) -- node[above, pos=0.5]{-1} (y);
			%		\draw[connection] (x2) -- node[above, pos=0.65]{1} (y);
			%		
			%	\end{tikzpicture}
		%
		\caption{A small NN with two input neurons $\bm x_1$ and $\bm x_2$, one
			hidden neuron labeled with the shape of the ReLU function, and one output neuron $\bm y$. The weights are given by the labels next to the arcs and all nodes have bias $0$. The NN has depth $2$, width $1$, and size~$1$. It
			computes the function $\bm x \mapsto \bm y= \bm x_2-\max\set{0,\bm x_2-\bm x_1}=-\max\set{-\bm x_2,-\bm x_1}
			=\min\set{\bm x_1,\bm x_2}$.}
		\label{Fig:Min2Num}
	\end{figure}
	
	\section{Max-Affine Arithmetic Programs}\label{Sec:MAAP}
	
	One way of specifying an NN is by explicitly writing down the network architecture, weights, and biases, that is, the affine transformations of each layer. However, for NNs that mimic the execution of an algorithm this is very unhandy and not well readable. For the purpose of an easier notation we introduce a pseudo-code language, called \emph{Max-Affine Arithmetic Programs (MAAPs)}, and prove that it is essentially equivalent to NNs.
	
	MAAPs perform arithmetic operations on real-valued \emph{variables}. In addition, there are natural- or real-valued \emph{constants} and different kinds of \emph{instructions}. In order to distinguish constants from variables the latter will be denoted by bold symbols. Each MAAP consists of a fixed number of input and output variables, as well as a sequence of (possibly nested) instructions.
	
	In order to describe an algorithm with an arbitrary number of input variables and to be able to measure asymptotic complexity, we specify a \emph{family} of MAAPs that is parametrized by a natural number $n$ that determines the number of input variables and is treated like a constant in each single MAAP of the family. A MAAP family then corresponds to a family of NNs. A similar concept is known from circuit complexity, where Boolean circuit families are used to measure complexity; cf.~\cite{arora2009computational}.
	
	MAAPs consist of the following types of instructions:
	
	\begin{enumerate}
		\item Assignment: this instruction assigns an expression to an old or new variable. The only two types of allowed expressions are affine combinations or maxima of affine combinations of variables:
		$b + \sum_{j} c_j\bm v_j$ and $\max \left\{\,b^{(i)} + \sum_{j} c_j^{(i)} \bm v_j^{(i)} \mathrel{}\middle|\mathrel{} i = 1, \dots, n\,\right\}$,
		where $n \in \N$, $b, b^{(i)}, c_j, c_j^{(i)} \in \R$ are constants and $\bm v_j, \bm v_j^{(i)} \in \R$ are variables.
		Without loss of generality minima are also allowed.
		\item \textbf{Do-Parallel}: this instruction contains a constant number of blocks of instruction sequences, each separated by an \textbf{and}. These blocks must be executable in parallel, meaning that each variable that is assigned in one block cannot appear in any other block.\footnote{Local variables in different blocks may have the same name if their scope is limited to their block.}
		%  For better readability we often omit the indent and the \textbf{Do-Parallel} key word. Instead an \textbf{and} keyword at the end of the line means that this instruction is done in parallel to the next instruction.
		\item \textbf{For-Do} loop: this is a standard for-loop with a \emph{constant number of iterations} that are executed sequentially.
		\item \textbf{For-Do-Parallel} loop: this is a for-loop with a \emph{constant number of iterations} in which the iterations are executed in parallel. Therefore, variables assigned in one iteration cannot be used in any other iteration.\footnote{Again, local variables within different iterations may have the same name if their scope is limited to their iteration.} 
	\end{enumerate}
	
	\Cref{Alg:Example} shows an example MAAP to illustrate the possible instructions.
	\begin{algorithm}[t]
		\SetKwFunction{Instructions}{Instructions}
		\SetKwInput{KwInput}{Input}
		\SetKwFor{ForParallel}{for each}{do parallel}{end}
		\SetKwIF{DoParallel}{AndIf}{And}{do parallel}{}{and}{and}{end}
		
		
		\caption{Instructions.\label{Alg:Example}}
		\KwInput{Input variables $\bm v_1, \bm v_2, \dots, \bm v_n$.}\vspace{0.5em}
		\tcp{Assignments and Expressions:}
		$\bm x_1 \leftarrow 4 +\sum_{i = 1}^n (-1)^i \cdot \bm v_i$ \\
		$\bm x_2 \leftarrow \max\set{3 \cdot \bm v_1, -1.5 \cdot \bm v_n, \bm x_1, 5}$
		\vspace{0.5em}
		
		\tcp{For-Do loop:} 
		\For{$k = 1, \dots, n-1$}{
			$\bm v_{k+1} \leftarrow \bm v_k + \bm v_{k+1}$
		} \vspace{0.5em}
		
		\tcp{Do-Parallel:}       
		\uDoParallel{}{
			$\bm y_1 \leftarrow \max \set{\bm x_1, \bm x_2}$}
		\uAnd{
			$\bm y_2 \leftarrow 7$
		} 
		\And{
			$\bm y_3 \leftarrow \sum_{i = 1}^n \bm v_i$
		} \vspace{0.5em}
		
		\tcp{For-Do-Parallel loop:} 
		\ForParallel{$k = 4, \dots, n$}{
			$\bm y_{k} \leftarrow \bm v_{k-1} - \bm v_{k}$\\
			$\bm y_k \leftarrow \max\set{\bm y_k, 0}$
		} \vspace{0.5em}
		
		\Return{$(\bm y_1, \bm y_2, \dots, \bm y_n)$}
	\end{algorithm}
	
	Note that we do not allow any \textbf{if}-statements or other branching operations. In other words, the number of executed instructions of an algorithm is always the same independent of the input variables.
	
	%Before we prove that a function can be computed by a MAAP if and only if it can be computed by an NN, that is, if and only if it is a CPWL function, we introduce three complexity measures $d(A)$, $w(A)$, and $s(A)$ for a MAAP~$A$. We will then see that they yield a direct correspondence to depth, width, and size of a corresponding NN. 
	In order to connect MAAPs with NNs, we introduce three complexity measures $d(A)$, $w(A)$, and $s(A)$ for a MAAP~$A$. We will then see that they yield a direct correspondence to depth, width, and size of a corresponding NN.
	
	For these complexity measures for MAAPs, assignments of affine transformations come ``for free'' since in an NN this can be realized ``between layers''. This is a major difference to other (parallel) models of computation, e.g., the parallel random access machine (PRAM)~\cite{greenlaw1995limits}. Apart from that, the complexity measures are recursively defined as follows.
	
	\begin{itemize}
		\item For an assignment $A$ with a maximum or minimum expression of $k\geq 2$ terms we have $d(A) \coloneqq \lceil \log_2 k \rceil$, $w(A)\coloneqq2k$, and $s(A)\coloneqq4k$.
		\item For a sequence $A$ of instruction blocks $B_1, B_2, \dots, B_k$ we have $d(A)\coloneqq\sum_{i=1}^k d(B_i)$, $w(A)\coloneqq\max_{i=1}^k w(B_i)$, and $s(A)\coloneqq\sum_{i=1}^k s(B_i)$.
		\item For a \textbf{Do-Parallel} instruction $A$ with parallel blocks $B_1, B_2, \dots, B_k$ we have \mbox{$d(A)\coloneqq\max_{i=1}^k d(B_i)$}, $w(A)\coloneqq\sum_{i=1}^k w(B_i)$, and $s(A)\coloneqq\sum_{i=1}^k s(B_i)$.
		\item For a \textbf{For-Do} loop $A$ with $k$ iterations that executes block $B_i$ in iteration $i$ we have $d(A)\coloneqq\sum_{i=1}^k d(B_i)$, $w(A)\coloneqq\max_{i=1}^k w(B_i)$, and $s(A)\coloneqq\sum_{i=1}^k s(B_i)$.
		\item For a \textbf{For-Do-Parallel} loop with $k$ iterations that executes block $B_i$ in iteration~$i$ we have $d(A)\coloneqq\max_{i=1}^k d(B_i)$, $w(A)\coloneqq\sum_{i=1}^k w(B_i)$, and $s(A)\coloneqq\sum_{i=1}^k s(B_i)$.
	\end{itemize}
	
	The following proposition establishes the desired correspondence between the complexities of MAAPs and NNs.
	
	\begin{restatable}{proposition}{NNMAAP}\label{Prop:NN-MAAP}
		For a function $f\colon\R^n\to\R^m$ the following is true.
		\begin{enumerate}[(i)]
			\item If $f$ can be computed by a MAAP $A$, then it can also be computed by an NN with depth $d(A)+1$, width $w(A)$, and size $s(A)$.
			\item If $f$ can be computed by an NN with depth $d+1$, width $w$, and size $s$, then it can also be computed by a MAAP $A$ with $d(A)=d$, $w(A)=2w$, and $s(A)=4s$.
		\end{enumerate}
	\end{restatable}

		\begin{proof}\
		\begin{enumerate}[(i)]
			\item First note that we can assume without loss of generality that $A$ does not contain \textbf{For-Do} or \textbf{For-Do-Parallel} loops. Indeed, since only a constant number of iterations is allowed in both cases, we can write them as a sequence of blocks or a \textbf{Do-Parallel} instruction, respectively. Note that this also does not alter the complexity measures $d(A)$, $w(A)$, and $s(A)$ by their definition. Hence, suppose for the remainder of the proof that $A$ consists only of assignments and (possibly nested) \textbf{Do-Parallel} instructions.
			
			The statement is proven by an induction on the number of lines of $A$. For the induction base suppose $A$ consists of a single assignment. If this is an affine expression, then an NN without hidden units (and hence with depth 1, width 0, and size 0) can compute $f$. If this is a maximum (or minimum) expression of $k$ terms, then, using the construction of~\cite[Lemma~D.3]{Arora:DNNwithReLU}, an NN with depth $\lceil \log_2 k \rceil + 1$, width $2k$, and size $4k$ can compute~$f$, which settles the induction base.
			
			For the induction step we consider two cases. If $A$ can be written as a sequence of two blocks~$B_1$ and~$B_2$, then, by induction, there are two NNs representing~$B_1$ and~$B_2$ with depth $d(B_i)+1$, width $w(B_i)$, and size $s(B_i)$ for $i=1,2$, respectively. An NN representing $A$ can be obtained by concatenating these two NNs, yielding an NN with depth $d(B_1)+d(B_2)+1=d(A)+1$, width $\max\set{w(B_1), w(B_2)} = w(A)$, and size~$s(B_1)+s(B_2)=s(A)$; cf.~\cite[Lemma~D.1]{Arora:DNNwithReLU}. Otherwise, $A$ consists of a unique outermost \textbf{Do-Parallel} instruction with blocks $B_1, B_2, \dots, B_k$. By induction, there are $k$ NNs representing $B_i$ with depth $d(B_i)+1$, width $w(B_i)$, and size~$s(B_i)$,~$i\in[k]$, respectively. An NN representing $A$ can be obtained by plugging all these NNs in parallel next to each other, resulting in an NN of depth $\max_{i=1}^k d(B_i)+1 = d(A) + 1$, width $\sum_{i=1}^k w(B_i) = w(A)$, and size $\sum_{i=1}^k s(B_i) = s(A)$. This completes the induction.
			
			\item Suppose the NN is given by a directed graph $G=(V,E)$ as described in \Cref{Sec:Prelim}. It is easy to verify that \Cref{Alg:NNtoMAAP} computes $f$ with the claimed complexity measures.
			\qedhere
		\end{enumerate}
	\end{proof}
	
	\begin{algorithm}[t]
		\SetKwInput{KwInput}{Input}
		\SetKwFor{ForParallel}{for each}{do parallel}{end}
		
		\caption{A generic MAAP to execute a given NN.\label{Alg:NNtoMAAP}}
		\KwInput{Input variables $\bm o(v)$ for $v\in V_0$.}\vspace{0.5em}
		
		\tcp{For each hidden layer:} 
		\For{$\ell = 1, \dots, d$}{
			\tcp{For each neuron in the layer:}
			\ForParallel{$v\in V_\ell$}{
				$\bm a(v) \leftarrow b_v + \sum_{u\colon uv\in E} w_{uv} \bm o(u)$\\
				$\bm o(v) \leftarrow \max\set{0, \bm a(v)}$
			}
		} \vspace{0.5em}
		
		\tcp{For each output neuron:} 
		\ForParallel{$v\in V_{d+1}$}{
			$\bm a(v) \leftarrow b_v + \sum_{u\colon uv\in E} w_{uv} \bm o(u)$
		} \vspace{0.5em}
		
		\Return{$(\bm a(v))_{v\in V_{d+1}}$.}
	\end{algorithm}
	
	\section{Exact Neural Networks for Minimum Spanning Trees}
	\label{Sec:CO}
	
	A spanning tree in a graph is a set of edges that is connected, spans all vertices, and does not contain any cycle. For given edge weights, the Minimum Spanning Tree Problem is to find a spanning tree with the least possible total edge weight.

	Most common algorithms, for example Kruskal's or Prim's algorithm, compare the edge weights and use conditional branchings to determine the order in which edges are added to the solution. Thus, they cannot be written as a MAAP or implemented as NN. Nevertheless, by ``tropicalizing'' a result by Fomin et al.~\cite{fomin2016subtraction} from arithmetic circuit complexity, we obtain the following proposition.
	
	\begin{restatable}{proposition}{propMST}\label{Prop:MST}
		For a fixed graph with $n$ vertices, there exists an NN of depth $\mathcal{O}(n\log n)$, width $\mathcal{O}(n^2)$, and size $\mathcal{O}(n^3)$ that correctly maps a vector of edge weights to the value of a minimum spanning tree.
	\end{restatable}
	
	To be more precise, Fomin et al.~\cite{fomin2016subtraction} provide a construction of a polynomial-size \emph{subtraction-free arithmetic circuit} (with standard addition, multiplication, and division, but without subtractions) to compute the so-called \emph{spanning tree polynomial} of a graph $(V,E)$, that is, the polynomial
	\[
	\sum_{T\text{ spanning tree}}\ \prod_{e\in T} x_e
	\]
	defined over $\abs{E}$ many variables $x_e$ associated with the edges of the graph.
	
	Tropicalizing this polynomial (to min-plus algebra) results precisely in the tropical polynomial mapping edge weights to the value of a minimum spanning tree:
	\[
	\min_{T\text{ spanning tree}}\ \sum_{e\in T} x_e.
	\]
	
	In the same way, one can tropicalize the arithmetic circuit provided by Fomin et al.~\cite{fomin2016subtraction}. In fact, every sum gate is just replaced with the small NN from \Cref{Fig:Min2Num} computing the minimum of its inputs, every product with a summation, and every division with a subtraction (realized using negative weights). That way, we obtain a polynomial-size NN to compute the value of a minimum spanning tree from any given edge weights. Note that it is crucial that the circuit is \emph{subtraction-free} because there is no inverse with respect to tropical addition.
	
	While this argument is already sufficient to justify the existence of polynomial-size NNs to compute the value of a minimum spanning tree, we now give a completely combinatorial proof of \Cref{Prop:MST}. The resulting NN, however, is basically equivalent to the NN obtained by tropicalizing the circuit of Fomin et al.~\cite{fomin2016subtraction}.
	
	The rough idea of the proof is to use a recursive MAAP, which then translates to an NN of the required size. In each step of the recursion, one node of the graph is deleted and all remaining edge weights are updated in such a way that the objective value of the minimum spanning tree problem in the original graph can be calculated from the objective value in the smaller graph. This idea is basically the translation of the so-called star-mesh transformations used by Fomin~et~al.~\cite{fomin2016subtraction} into the combinatorial world.
	
	Without loss of generality, we restrict ourselves to complete graphs. Edges missing in the actual input graph can be represented with large weights such that they will never be included in a minimum spanning tree.
	
	For $n=2$ vertices, the MAAP simply returns the weight of the only edge of the graph. For $n\geq 3$, our MAAP is given in \Cref{Alg:MST}.
	
	\begin{algorithm}[t]
		\SetKwFor{ForParallel}{for each}{do parallel}{end}
		\SetKwInput{KwInput}{Input}
		\caption{$\MST_n$: Compute the value of a minimum spanning tree for the complete graph on $n\geq3$ vertices.} \label{Alg:MST}
		\KwInput{Edge weights $(\bm x_{ij})_{1\leq i < j \leq n}$.}\vspace{0.5em}
		$\bm y_n \leftarrow \min_{i\in[n-1]} \bm x_{in}$\label{Line:compy}\\
		\ForParallel{$1\leq i<j\leq n-1$}{
			$\bm x'_{ij} \leftarrow \min\set{\bm x_{ij},\ \bm x_{in} + \bm x_{jn} - \bm y_n}$ \label{Line:xprime}
		} \vspace{0.5em} 
		\Return{$\bm y_n + \MST_{n-1}\left((\bm x'_{ij})_{1\leq i < j \leq n-1}\right)$}
	\end{algorithm}
	
	Next, we prove the correctness of this MAAP, which allows us to prove \Cref{Prop:MST} afterwards.
	
	\begin{proposition}\label{Prop:MSTMAAP}
		\Cref{Alg:MST} correctly computes the value of a minimum spanning tree in the complete graph on $n$ vertices.
	\end{proposition}
	\begin{proof}
		We use induction on $n$. The trivial case $n=2$ settles the induction start. Now suppose that the subroutine $\MST_{n-1}$ correctly computes the value of an MST for $n-1$ vertices. We need to show that the returned value $\bm y_n + \MST_{n-1}\left((\bm x'_{ij})_{1\leq i < j \leq n-1}\right)$ is indeed the MST value for $n$ vertices.
		
		First, we show that the value computed by \Cref{Alg:MST} is not larger than the correct objective value. For this purpose, let $T$ be the set of edges corresponding to an MST of $G$. By potential relabeling of the vertices, assume that $\bm y_n=\bm x_{1n}$. Note that we may assume without loss of generality that $v_1v_n\in T$: if this is not the case, adding it to $T$ creates a cycle in $T$ involving a second neighbor $v_i\neq v_1$ of $v_n$. Removing $v_iv_n$ from $T$ results again in a spanning tree with total weight at most the original weight.
		
		We construct a spanning tree $T'$ of the subgraph spanned by the first $n-1$ vertices as follows: $T'$ contains all edges of $T$ that are not incident with $v_n$. Additionally, for each~$v_iv_n\in T$, except for $v_1v_n$, we add the edge $v_1v_i$ to $T'$. It is immediate to verify that this construction results in fact in a spanning tree.
		We then obtain
		
		\begin{align*}
			\sum_{v_iv_j\in T} \bm x_{ij}
			&= \bm x_{1n} + \sum_{v_iv_n\in T,\ i>1} \bm x_{in} + \sum_{v_iv_j\in T,\ i,j<n} \bm x_{ij}\\
			&= \bm y_n + \sum_{v_iv_n\in T,\ i>1} (\bm x_{in} + \bm x_{1n} - \bm y_n) + \sum_{v_iv_j\in T,\ i,j<n} \bm x_{ij}\\
			&\geq \bm y_n + \sum_{v_iv_n\in T,\ i>1} \bm x'_{1i} + \sum_{v_iv_j\in T,\ i,j<n} \bm x'_{ij}\\
			&= \bm y_n + \sum_{v_iv_j\in T'} \bm x'_{ij}\\
			&\geq \bm y_n + \MST_{n-1}\left((\bm x'_{ij})_{1\leq i < j \leq n-1}\right).
		\end{align*}
		
		Here, the first inequality follows by the way how the values of $\bm x'$ are defined in line~\ref{Line:xprime} and the second inequality follows since $T'$ is a spanning tree of the first $n-1$ vertices and, by induction, the MAAP is correct for up to $n-1$ vertices. This completes the first proof that the MAAP does not overestimate the objective value. 
		
		In order to show that the MAAP does not underestimate the true objective value, let~$T'$ be the set of edges of a minimum spanning tree of the first $n-1$ vertices. Let~\mbox{$E^*\subseteq T'$} be the subset of edges~$v_iv_j$, with \mbox{$1\leq i<j \leq n-1$}, in $T'$ that satisfy \mbox{$\bm x'_{ij}=\bm x_{in} + \bm x_{jn}-\bm y_n$}. Note that, in particular, we have $\bm x'_{ij}=\bm x_{ij}$ for all $v_iv_j\in T'\setminus E^*$, which will become important later. We show that we may assume without loss of generality that $E^*$ only contains edges incident with $v_1$. To do so, suppose there is an edge $v_iv_j\in E^*$ with $2\leq i<j \leq n-1$. Removing that edge from $T'$ disconnects exactly one of the two vertices $v_i$ and $v_j$ from $v_1$; say, it disconnects $v_j$. We then can add $v_1v_j$ to $T'$ and obtain another spanning tree in $G'$. Moreover, by the definition of the weights $\bm x'$ and the choice of~$v_1$, we obtain $\bm x'_{1j}\leq \bm x_{1n} + \bm x_{jn} - \bm y_n \leq \bm x_{in} + \bm x_{jn}-\bm y_n = \bm x'_{ij}$. Hence, the new spanning tree is still minimal. This procedure can be repeated until every edge in $E^*$ is incident with $v_1$.
		
		Now, we construct a spanning tree $T$ in $G$ from $T'$ as follows: $T$ contains all edges of~$T'\setminus E^*$. Additionally, for every $v_1v_i\in E^*$, we add the edge $v_iv_n$ to $T$. Finally, we also add~$v_1v_n$ to $T$. Again it is immediate to verify that this construction results in fact in a spanning tree, and we obtain
		
		\begin{align*}
			\sum_{v_iv_j\in T} \bm x_{ij}
			&= \bm x_{1n} + \sum_{v_1v_i\in E^*} \bm x_{in} + \sum_{v_iv_j\in T'\setminus E^*} \bm x_{ij}\\
			&= \bm y_n + \sum_{v_1v_i\in E^*} (\bm x_{in} + \bm x_{1n} - \bm y_n) + \sum_{v_iv_j\in T'\setminus E^*} \bm x_{ij}\\
			&= \bm y_n + \sum_{v_1v_i\in E^*} \bm x'_{1i} + \sum_{v_iv_j\in T'\setminus E^*} \bm x'_{ij}\\
			&= \bm y_n + \sum_{v_iv_j\in T'} \bm x'_{ij}\\
			&= \bm y_n + \MST_{n-1}\left((\bm x'_{ij})_{1\leq i < j \leq n-1}\right).
		\end{align*}
		
		This shows that the MAAP returns precisely the value of the spanning tree $T$. Hence, its output is at least as large as the value of an MST, completing the second direction.
	\end{proof}
	
	Finally, we prove complexity bounds for the MAAP, allowing us to bound the size of the corresponding NN.
	
	\propMST*
	\begin{proof}
		In \Cref{Prop:MSTMAAP}, we have seen that \Cref{Alg:MST} performs the required computation. We show that $d(\MST_n)=\mathcal{O}(n\log n)$, $w(\MST_n)=\mathcal{O}(n^2)$, and $s(\MST_n)=\mathcal{O}(n^3)$. Then, the claim follows by \Cref{Prop:NN-MAAP}.
		
		Concerning the complexity measure $d$, observe that in each recursion the bottleneck is to compute the minimum in line~\ref{Line:compy}. This is of logarithmic order. Since we have $n$ recursions, it follows that $d(\MST_n)=\mathcal{O}(n\log n)$.
		
		Concerning the complexity measure $w$, observe that the bottleneck is to compute the parallel \textbf{for} loop in line~\ref{Line:xprime}. This is of quadratic order, resulting in $w(\MST_n)=\mathcal{O}(n^2)$.
		
		Finally, concerning the complexity measure $s$, the bottleneck is also the parallel \textbf{for} loop in line~\ref{Line:xprime}. Again, this is of quadratic order and since we have $n$ recursions, we arrive at~$s(\MST_n)=\mathcal{O}(n^3)$.
	\end{proof}
	
	
	
	\section{Exact Neural Networks for Maximum Flows}\label{Sec:MaxFlow}
	
	In this section we show that, given a fixed directed graph with $n$ nodes and $m$ edges, there is a polynomial-size NN computing a function of type~\mbox{$\R^{m}\to\R^{m}$} that maps arc capacities to a corresponding maximum flow. We achieve this by providing a MAAP family of polynomial complexity and using \Cref{Prop:NN-MAAP}.

	\mypara{The Maximum Flow Problem}
	
	Let $G = (V, E)$ be a directed graph with a finite node set $V =
	\set{v_1, \dots, v_n}$, $n \in \N$, containing a source $s = v_1$, a sink $t = v_n$, and an arc set $E \subseteq
	V^2 \setminus \set{vv|v \in V}$ in which each arc~$e \in E$ is equipped with a capacity $\nu_e \geq 0$. We write $m=\abs{E}$ for the number of arcs,~$\delta^+_v$~and~$\delta^-_v$ for the outgoing and incoming arcs of node $v$, as well as $N^+_v$ and~$N^-_v$ for the successor and predecessor nodes of $v$ in $G$, respectively. The distance~$\dist_{G}(v,w)$ denotes the minimum number of arcs on any path from $v$ to $w$ in $G$.
	
	The \emph{Maximum Flow Problem} consists of finding an $s$-$t$-flow $(y_e)_{e \in
		E}$ satisfying $0 \leq y_e \leq \nu_e$ and~\mbox{$\sum_{e \in \delta^-_v} y_e = \sum_{e \in \delta^+_v} y_e$} for all $v \in V
	\setminus \set{s,t}$ such that the \emph{flow value} $\sum_{e \in \delta^+_s} y_e - \sum_{e \in \delta^-_s} y_e$ is maximal.
	
	For the sake of an easier notation we assume for each arc $e=uv\in E$ that its reverse arc~$vu$ is also contained in $E$. This is without loss of generality because we can use capacity~$\nu_e=0$ for arcs that are not part of the original set $E$.
	%
	In order to avoid redundancy we represent flow only in one arc direction. More precisely, with $\vec E=\set{v_iv_j\in E | i<j}$ being the set of \emph{forward arcs}, we denote a flow by $(y_e)_{e \in \vec E}$. The capacity constraints therefore state that $-\nu_{vu} \leq y_{uv} \leq \nu_{uv}$. Hence, a negative flow value on a forward arc $uv \in \vec E$ denotes a positive flow on the corresponding \emph{backward arc} $vu$.
	
	A crucial construction for maximum flow algorithms is the \emph{residual network}. For a given $s$-$t$-flow $(y_{e})_{e \in \vec E}$, the \emph{residual capacities} are defined as follows. For an arc $uv \in \vec E$ the \emph{residual forward capacity} is given by $c_{uv} \coloneqq \nu_{uv} - y_{uv}$ and the \emph{residual backward capacity} by~$c_{vu} \coloneqq \nu_{vu} + y_{uv}$. The \emph{residual network} consists of all directed arcs with positive residual capacity. Hence, it is given by $G^*=(V, E^*)$ with $E^* \coloneqq \set{e \in E | c_e > 0}$.
	
	\mypara{A MAAP to solve the Maximum Flow Problem}
	
	Our MAAP, given in \Cref{Alg:max_flow}, uses ideas that are similar to the algorithms by Dinic~\cite{dinic1970algorithm} and by Edmonds and Karp~\cite{edmonds1972theoretical}. We start with the zero flow and consider the residual network $G^*$. In each iteration we augment flow along arcs of the residual network that lie on a shortest path from the source $s$ to the sink $t$ until they become disconnected in $G^*$. However, in contrast to the classical algorithms, finding an appropriate augmenting flow is much more technically involved, due to the limited set of operations allowed. This is accomplished by the \texttt{FindAugmentingFlow$_k$} subroutine defined in \Cref{Alg:augmenting_flow}. Its key feature is to return a flow that has positive flow values only on arcs that lie on a path of length $k$ from $s$ to $t$ in the current residual network. Moreover, if such a path still exists, at least one arc of the residual network is saturated by that flow. Of course, all this needs to happen without explicitly knowing the arcs contained in the current residual network since this would involve conditional branchings.
	
	\begin{algorithm}[t]
		%\SetKwFunction{MaxFlow}{MaxFlow}
		\SetKwFunction{FindAugmentingFlow}{FindAugmentingFlow$_k$}
		\SetKwFor{ForParallel}{for each}{do parallel}{end}
		\SetKwInput{KwInput}{Input}
		\caption{Compute a maximum flow for a fixed graph~\mbox{$G=(V,E)$}.} \label{Alg:max_flow}
		\KwInput{Capacities $(\bm \nu_e)_{e \in E}$.}\vspace{0.5em}
		\tcp{Initializing:}
		\ForParallel{$uv\in \vec E$}{\label{Line:startIni}
			$\bm x_{uv} \leftarrow 0$ %\textbf{ and} \hspace{9.5pt}%
			\hspace{9.5pt}\tcp{flow; negative value correspond to flow on $vu$} \label{Line:ini1}
			$\bm c_{uv} \leftarrow \bm \nu_{uv}$ %\textbf{ and}%
			\tcp{residual forward capacities}
			$\bm c_{vu} \leftarrow \bm \nu_{vu}$ %\hspace{22pt} 
			\tcp{residual backward capacities}\label{Line:endIni}
		} \vspace{0.5em}           
		\tcp{Main part:} 
		\For{$k = 1, \dots, n-1$}{\label{Line:outerFor}
			\For{$i = 1, \dots, m$}{\label{Line:innerFor}
				$(\bm y_e)_{e \in \vec E} \leftarrow \FindAugmentingFlow((\bm c_e)_{e \in E})$ \label{Line:sub}\\
				\tcc{Returns an augmenting flow (respecting the residual capacities) that only uses paths of length exactly $k$ and saturates at least one arc.}
				\tcp{Augmenting:}
				\ForParallel{$uv \in \vec E$}{\label{Line:startUpdate}
					$\bm x_{uv} \leftarrow \bm x_{uv} + \bm y_{uv}$% \textbf{ and} 
					\label{Line:beginUpdate}\\
					$\bm c_{uv} \leftarrow \bm c_{uv} - \bm y_{uv}$% \textbf{ and} 
					\\
					$\bm c_{vu} \leftarrow \bm c_{vu} + \bm y_{uv}$\label{Line:endUpdate}
				}
			}
		}
		\Return{$(\bm x_{e})_{e \in \vec E}$}
	\end{algorithm}

	Note that there are several instances in the MAAP where instructions could be run in parallel (e.g., line~\ref{Line:ini1} to \ref{Line:endIni} or line~\ref{Line:beginUpdate} to \ref{Line:endUpdate}). As we only strive for asymptotic complexity bounds and since running these instructions in parallel would only change constant factors, we omit the \textbf{Do-Parallel} statements for the sake of readability.	

	\mypara{The FindAugmentingFlow Subroutine}
	
	The key component of the Maximum-Flow-MAAP is the subroutine, given by \Cref{Alg:augmenting_flow}, that returns an augmenting flow which uses only paths of at most a fixed length $k$ and, in addition, saturates at least one arc of the residual network.
	
	\begin{algorithm}[tp]
		\SetKwInput{KwInput}{Input}
		\SetKwFor{ForParallel}{for each}{do parallel}{end}
		\caption{FindAugmentingFlow$_k$ for a fixed graph $G=(V,E)$ and a fixed length $k$.} \label{Alg:augmenting_flow}
		\KwInput{Residual capacities $(\bm c_e)_{e \in E}$.}
		\vspace{0.5em}
		
		\tcp{Initializing:}
		\ForParallel{$vw \in \vec E$}{\label{line:start}
			$\bm z_{vw} \leftarrow 0$ % \textbf{ and}  
			\quad \tcp{flow in residual network}
			$\bm z_{wv} \leftarrow 0$
		}
		\ForParallel{$(i, v) \in [k] \times  (V \setminus \set{t})$}{
			$\bm Y_v^i \leftarrow 0$ % \textbf{ and} \quad 
			\hspace{3pt}\tcp{excessive flow at $v$ in iteration $i$ (from $k$ to $1$)}
			$\bm a_{i,v}\leftarrow 0$ \label{line:init_fattest}\tcp{initialize fattest path values}
		} 
		\vspace{0.5em}
		
		\tcp{Determining the fattest path values:} 
		\ForParallel{$v\in N^-_t$}{ \label{line:begin_fattest_path}
			$\bm a_{1, v} \leftarrow \bm c_{vt}$\label{line:fattest_t}
		}
		\For{$i = 2, 3, \dots, k$}{\label{line:begin_fattest_for}
			\ForParallel{$v \in V \setminus\set{t}$}{
				$\bm a_{i, v} \leftarrow \max_{w \in N^+_v \setminus\set{t}}\min\set{\bm a_{i-1,w}, \bm c_{vw}}$ \label{line:end_fattest_path}
			}            
		}
		\vspace{0.5em}
		
		\tcp{Pushing flow of value $\bm a_{k, s}$ from $s$ to $t$:}
		$\bm Y_s^k \leftarrow \bm a_{k, s}$ \tcp{excessive flow at $s$} \label{line:begin_pushing}
		\For{$i = k, k-1, \dots, 2$}{\label{line:outer_for}
			\For{$v \in V\setminus\set{t}$ in index order}{
				\For{$w\in N^+_v\setminus\set{t}$ in index order}{
					\tcp{Push flow out of $v$ and into $w$:}
					$\bm f \leftarrow \min\set{\bm Y_v^i, \bm c_{vw}, \bm a_{i-1, w} - \bm Y_w^{i-1}}$ \tcp{value we can push over $vw$ such that this flow can still arrive at $t$}\label{line:find_push_amount}
					$\bm z_{vw} \leftarrow \bm z_{vw} + \bm f$ %\textbf{ and}
					\label{line:increase_of_z}\\
					$\bm Y_v^i \leftarrow \bm Y_v^i - \bm f$ %\textbf{ and}
					\\
					$\bm Y_w^{i-1} \leftarrow \bm Y_w^{i-1} + \bm f\label{line:increase_of_Y}$
				}
			}
		}
		\ForParallel{$v\in N^-_t$}{\label{line:start_t_pushing}
			\tcp{Push flow out of $v$ and into $t$:}
			$\bm z_{vt} \leftarrow \bm Y_v^1$\\
			$\bm Y_v^1 \leftarrow 0$ 
		} \label{line:end_pushing}
		\vspace{0.5em}
		
		\tcp{Clean-up by bounding:}
		\For{$i = 2, 3, \dots, k - 1$}{ \label{line:begin_clean_up}
			\For{$w \in V \setminus \set{t}$ in reverse index order}{    
				\For{$v \in N^-_w \setminus \set{t}$ in reverse index order}{
					$\bm b \leftarrow \min\set{\bm Y_w^i, \bm z_{vw}}$ \tcp{value we can push backwards along $vw$}
					$\bm z_{vw} \leftarrow \bm z_{vw} - \bm b$ %\textbf{ and}
					\\
					$\bm Y_w^i \leftarrow \bm Y_w^i - \bm b$ %\textbf{ and}
					\\
					$\bm Y_v^{i+1} \leftarrow \bm Y_v^{i+1} + \bm b$    \label{line:increase_of_Y_by_b} \label{line:end_clean_up}       
				}
			}
		}
		\vspace{0.5em}
		\ForParallel{$uv \in \vec E$}{\label{line:start_y}
			$\bm y_{vw} \leftarrow \bm z_{vw} - \bm z_{wv}$ \label{line:end}
		}
		
		\Return{$(\bm y_e)_{e \in \vec E}$}
	\end{algorithm}

	The first step of this subroutine is to determine for each node $v$ and each $i \in \N_0$ the maximal flow value $\bm a_{i, v}$ that can be
	sent from $v$ to $t$ on a single path of length exactly~$i$ in the residual network. We call such a path a
	\emph{fattest path} of length $i$ from $v$ to $t$.
	The value~$\bm a_{k,s}$ of a fattest path from $s$ to $t$ of length $k$ is of particular interest as the algorithm greedily pushes
	this value from $s$ to $t$ within $k$ iterations. This means that at each node $v$ flow is pushed into arcs $vw \in
	\delta_v^+$ in the order given by the node indices of the successors $w$. Here, the flow value that is pushed into a node $w$ should
	not exceed $\bm a_{i,w}$, where $i$ is the current iteration. This way, all flow that is pushed into a node could in
	principle reach $t$ within $i$ steps.
	However, it can happen that not all the flow that is pushed into a node $v$ can be pushed out of $v$ since the next nodes might be loaded with flow already. Therefore, the inflow at some nodes might be larger than the outflow after the pushing procedure; see \Cref{fig:find_augmenting_flow} for an example where this happens.
	
	In order to restore flow conservation the algorithm performs a clean-up in which it iterates over the nodes in reverse order and reduces
	the incoming flow by the remaining excess flow at node $v$. Hence, we obtain a feasible augmenting $s$-$t$-flow in the end. 
	%In \Cref{Alg:augmenting_flow}, we give a formal description of this subroutine and an example is given in \Cref{fig:find_augmenting_flow}.
	
	\begin{figure}[tp]
		\centering \includegraphics[page=1]{find_augmenting_flow}\\[2em]
		\includegraphics[page=2,width=0.49\linewidth]{find_augmenting_flow}\hfill
		\includegraphics[page=3,width=0.49\linewidth]{find_augmenting_flow}\\[1em]
		\includegraphics[page=4,width=0.49\linewidth]{find_augmenting_flow}\hfill
		\includegraphics[page=5,width=0.49\linewidth]{find_augmenting_flow}\\[2em]
		\includegraphics[page=6,width=0.49\linewidth]{find_augmenting_flow}\hfill
		\includegraphics[page=7,width=0.49\linewidth]{find_augmenting_flow}
		\caption{Example of the FindAugmentingFlow Subroutine. The network on the top depicts the residual capacity bounds $[-\bm c_{wv}, \bm c_{vw}]$ for all arcs $vw \in \vec E$ and the node labels $\bm a_{i,v}$ for the fattest path from $v$ to $t$ within $i$ steps. The four figures in the middle show the states of the flow $\bm z_{vw}$ and the excessive flow $\bm Y_v^i$ at the end of a push iteration. 
			The bottom two figures depict the states after the clean-up iterations. All values that are not displayed are zero.
			Observe that the result is an $s$-$t$-flow that is feasible with respect to the residual capacities, uses only paths of length $k=4$, and saturates the arc $v_6t$.}
		\label{fig:find_augmenting_flow}
	\end{figure}
	
	The following theorem states that the subroutine indeed computes an augmenting flow fulfilling all required properties.
	\begin{samepage}
		\begin{restatable}{theorem}{thmaugmentingflow} \label{Thm:augmenting_flow}
		Let $(\bm c_e)_{e \in E}$ be residual capacities such that the distance between $s$ and $t$ in the residual network $G^*=(V,
		E^*)$ is at least $k$. Then the MAAP given in \Cref{Alg:augmenting_flow} returns an $s$-$t$-flow $\bm y = (\bm y_{uv})_{{uv} \in \vec E}$ with $-\bm c_{vu} \leq \bm y_{uv} \leq \bm c_{uv}$ such that we have positive flow only on arcs that lie on an $s$-$t$-path of length exactly $k$ in $G^*$. If the distance of $s$ and $t$ in the residual network is exactly $k$, then $\bm y$
		has a strictly positive flow value and there exists at least one saturated arc, i.e., one arc $e \in E^*$ with
		$\bm y_{e} = \bm c_{e}$.
	\end{restatable}
	\end{samepage}
	The proof idea is as follows. For the flow conservation we first show that the $\bm Y_u$ variables do indeed track the excessive flow at node $u$. More precisely, we show that after the pushing procedure we have
	$\sum_{e \in \delta_u^-} \bm z_e - \sum_{e \in \delta_u^+} \bm z_e = \bm Y_u^{k - \dist_{G^*}(s,u)}$
	and $\bm Y_u^i = 0$ for all~\mbox{$i \neq k - \dist_{G^*}(s,u)$}. The clean-up reduces the excessive flow to zero, i.e., in the end it holds that $\bm Y_u^i = 0$ for all $i$ including $i = k - \dist_{G^*}(s,u)$.
	
	In order to show that at least one residual arc is saturated we consider a node $v^*$ that has positive excess flow
	after the pushing phase. Among these, $v^*$ is chosen as one of the closest nodes to~$t$ in~$G^*$. From all shortest $v^*$-$t$-paths in~$G^*$ we pick the path~$P$
	that has lexicographically the smallest string of node indices. As the fattest path from $v^*$ to $t$ has at least the
	residual capacity of $P$ (given by the minimal residual capacity of all arcs along $P$), the pushing procedure has pushed
	at least this value along $P$. Hence, the arc on $P$ with minimal capacity has to be saturated. It is then easy to show
	that the clean-up does not reduce the value along $P$.
	
	\begin{proof}[Proof of \Cref{Thm:augmenting_flow}]
		It is easy to check that lines~\ref{line:begin_fattest_path} to \ref{line:end_fattest_path} from the algorithm do
		indeed compute the maximal flow value $\bm a_{i,v}$ that can be send from $v$ to $t$ along a single path (which we call
		the \emph{fattest path}) of length exactly $i$.
		
		In the following we show that in line \ref{line:end} the arc vector $\bm z = (\bm z_e)_{e \in E}$ forms an $s$-$t$-flow in the residual network $G^* = (V,E^*)$ that satisfies $0 \leq \bm z_e \leq \bm c_e$.
		For this, recall that $\dist_{G^*}(s,u)$ denotes the distance from $s$ to $u$ in the residual network.
		
		In order to prove flow conservation of $\bm z$ at all vertices except for $s$ and $t$, we fix some node $u \in V \setminus \set{t}$ and show that
		\begin{equation} \label{eq:Y_is_excess_flow}  
			\bm Y_u^{j} = \begin{cases}
				\sum_{e \in \delta_u^-} \bm z_e - \sum_{e \in \delta_u^+} \bm z_e & \text{if } j = k - \dist_{G^*}(s,u),\\
				0 & \text{otherwise,}
			\end{cases}
		\end{equation}
		holds throughout the execution of the subroutine.
		
		\begin{claim} \label{claim:Y_is_excess}
			\Cref{eq:Y_is_excess_flow} holds after the pushing procedure (lines \ref{line:begin_pushing} to \ref{line:end_pushing}).
		\end{claim}
		
		\begin{subproof}[Proof of \Cref{claim:Y_is_excess}]
			For $j < k-\dist_{G^*}(s,u)$, there does not exist any $u$-$t$-path of length~$j$ (since $j + \dist_{G^*}(s,u) < k \leq \dist_{G^*}(s,t)$). Hence, the fattest
			$u$-$t$-path of length exactly~$j$ has capacity~$\bm a_{j,u}=0$.
			In any iteration that might increase $\bm Y_u^j$, that is, for $i=j+1$, $v \in V \setminus \set{u, t}$, and $w=u$, we have $\bm f = 0$. This implies that $\bm Y_u^j$ remains $0$.
			
			For $j > k - \dist_{G^*}(s,u)$, there is no $s$-$u$-path of length $k - j$ in the residual graph as $k - j < \dist_{G^*}(s,u)$.
			The statement $\bm Y_u^j = 0$ then follows by an induction on $\dist_{G^*}(s,u)$, as we show in the following:
			
			For the base case of $\dist_{G^*}(s,u) = 1$ we have that $\bm Y_u^j = \bm Y_u^k = 0$ since for $u \neq s$ it holds that $\bm Y_u^k$	stays $0$ during the whole algorithm.
			
			The induction step follows because for iteration $i = j+1$, $v \in V \setminus \set{u, t}$ and $w =
			u$ it holds that either $\bm c_{vw} = 0$ (i.e., arc $vw$ is not part of the residual network) or $\bm Y_v^i = 0$
			(inductively as $\dist_{G^*}(s,v) \geq \dist_{G^*}(s,u)-1$ and $i = j+1 > k - \dist_{G^*}(s,u) + 1 \geq k - \dist_{G^*}(s,v)$). Either way we have $\bm f = 0$ implying that
			$\bm Y_u^j = \bm Y_w^{i - 1}$ stays $0$.
			
			In conclusion, we obtain that $\bm Y_u^j$ can only be non-zero for $j = k-\dist_{G^*}(s,u)$. In each iteration
			with $i = k-\dist_{G^*}(s,u) +1$ and $w = u$ we add $\bm f$ to the flow value $\bm z_{vu}$ and the same to $\bm
			Y_u^{k-\dist_{G^*}(s,u)}$ and in each iteration with $i = k-\dist_{G^*}(s,u)$ and $v = u$ we add~$\bm f$ to the flow
			value $\bm z_{uw}$ and subtract $\bm f$ from $\bm Y_u^{k-\dist_{G^*}(s,u)}$. Hence, $\bm Y_u^{k-\dist_{G^*}(s,u)}$
			denotes exactly the excessive flow after the pushing procedure as stated in \eqref{eq:Y_is_excess_flow}.
		\end{subproof}
		
		This claim already shows that $\bm z_e$ can only be positive if $e$ lies on an $s$-$t$-path of length exactly $k$, which is a
		shortest path in the residual network. To see this, let $vw$ be an arc that is not on such a path. In line~\ref{line:find_push_amount}, it either holds that $\bm Y_v^i = 0$ (if $i \neq k - \dist_{G^*}(s,u)$) or $\bm a_{i - 1, w} = 0$ because for $i = k -
		\dist_{G^*}(s,u)$ there is no $w$-$t$-path of length $i-1 = k - \dist_{G^*}(s,u) - 1$ (since otherwise  $vw$ would
		lie on an $s$-$t$-path of length $k$). Thus, $\bm z_{vw}$ will never be increased. As the clean-up only reduces the flow values, $\bm z_{vw}$ will still be $0$ at the end (line \ref{line:end}).
		
		\begin{claim} \label{claim:Y_is_excess_after_clean_up}		\Cref{eq:Y_is_excess_flow} holds in each iteration of the clean-up (lines \ref{line:begin_clean_up} to \ref{line:end_clean_up}).
		\end{claim}
		\begin{subproof}[Proof of \Cref{claim:Y_is_excess_after_clean_up}]
			First, we show that $\bm Y_u^j$ stays $0$ for $j \neq k - \dist_{G^*}(s,u)$ by induction over $j = 2,3, \dots, k$. The
			base case follows immediately as we only subtract $\bm b \geq 0$ from $\bm Y_u^2$. For the induction step we have to
			show that $\bm b = 0$ whenever we add $\bm b$ to $\bm Y_u^j$ in line \ref{line:increase_of_Y_by_b}. In all iterations
			with $i = j-1$ and $v = u$ we either have $\bm z_{uw} = 0$ or $\bm Y_w^i = 0$. The reason for this is that $\bm z_{uw}
			> 0$ implies that $uw$ lies on a shortest $s$-$t$-path, which means that $\dist_{G^*}(s,w) = \dist_{G^*}(s,u) + 1$, and
			hence, $i = j-1 \neq k - \dist_{G^*}(s,u) - 1 = k - \dist_{G^*}(s,w)$. By induction this means that $\bm Y_w^i = 0$.
			Either way this implies $\bm b = 0$.
			
			\Cref{eq:Y_is_excess_flow} holds for $j = k - \dist_{G^*}(s,u)$ since for $e \in \delta_u^-$ the value $\bm b$ is only
			possibly positive for $i = k - \dist_{G^*}(s,u)$ and then it is subtracted from $\bm z_e$ as well as from $\bm Y_u^{k -\dist_{G^*}(s,u)}$. For $e
			\in \delta_u^+$, the value $\bm b$ can only be positive for $i = k - \dist_{G^*}(s,u) + 1$, and hence, $\bm b$ is subtracted from $\bm z_e$
			exactly when it is added to $\bm Y_u^{k - \dist_{G^*}(s,u)}$.
		\end{subproof}
		
		Next, we show that at the end of the subroutine it holds that $\bm Y_u^j = 0$ for all $j$, in particular also for $j = k -\dist_{G^*}(s,u)$. The only exception of this is
		$\bm Y_s^k$. To see this, first observe that during the clean-up, $\bm Y_u^{k - \dist_{G^*}(s,u)}$ is maximal after iteration $i = k
		- \dist_{G^*}(s,u)-1$ and does not increase anymore for $i \geq k - \dist_{G^*}(s,u)$. At the start of iteration $i = k - \dist_{G^*}(s,u)$ it holds due to \eqref{eq:Y_is_excess_flow} that
		\[\sum_{e \in \delta_u^-} \bm z_e \geq \bm Y_u^{k - \dist_{G^*}(s,u)}.\]
		Hence, for $i = k - \dist_{G^*}(s,u)$ and $w = u$ there is one iteration for all $e \in \delta_u^-$ and within this
		iteration $\bm Y_u^{k - \dist_{G^*}(s,u)}$ is reduced by $\bm z_e$ until $\bm Y_u^{k - \dist_{G^*}(s,u)} = 0$. This shows
		that after all iterations with $i = k - \dist_{G^*}(s,u)$ it holds that $\bm Y_u^{k - \dist_{G^*}(s,u)} = 0$. Together with~\eqref{eq:Y_is_excess_flow}, this immediately
		implies flow conservation of $(\bm z_e)_{e \in E}$.
		
		Finally, in order to show that $0 \leq \bm z_e \leq \bm c_e$, note that $\bm z_e$ is initialized with $0$ and it is only
		increased in line~\ref{line:increase_of_z} of the unique iteration with $vw = e$ and $i = k-\dist_{G^*}(s,v)$, as we have argued in
		the proof of \Cref{claim:Y_is_excess}. In this iteration we have that $\bm f \leq \bm c_e$, which immediately shows that
		$0 \leq \bm z_e \leq \bm c_e$.
		
		It only remains to show that at least one residual arc is saturated. To this end, suppose that the distance
		of $s$ and $t$ in $G^*$ is $k$, which means that there exists at least one $s$-$t$-path of length exactly $k$ with a strictly
		positive residual capacity on all arc along this path.
		
		Let us consider the set $\set{(v, i) |\bm Y^i_v  > 0 \text{ after the pushing procedure}}$. These are all nodes that need to be
		cleaned up in order to restore flow conservation, paired with their distance to $t$. Let $(v^*, i^*)$ be a tuple of this set such that $i^*$ is minimal.
		In other words, $v^*$ is a node that is closest to $t$ among theses nodes. 
		%This set should never be empty as $(s, \bm Y^k_s)$ is always included. But as we do not want to prove $\bm Y^k_s >0$, 
		We manually set $(v^*, i^*)$ to $(s, k)$ in the case that the set is empty.
		
		\begin{claim} \label{claim:sat_one} Some arc on a shortest path from $v^*$ to $t$ in $G^*$ is saturated by $\bm y_e$.
		\end{claim}
		\begin{subproof}[Proof of \Cref{claim:sat_one}]
			Among all these paths between $v^*$ and $t$ of length $i^*$ we consider the path $P$ which has lexicographically the
			smallest string of node indices. Let $\bm c_{\min}$ be the minimal residual capacity along this path $P$.
			
			For all nodes $v$ along $P$ (including $v^*$) we have $\bm a_{i,v} \geq \bm c_{\min}$, where $i$ is the distance from~$v$ to $t$ along $P$, since the fattest path from $v$ to $t$ has to have at least the residual capacity of $P$.
			
			After the pushing procedure it holds that $\bm z_e \geq \bm c_{\min}$ for all $e \in P$. This is true for the first arc
			on $P$, since we have excess flow at node $v^*$ remaining (after pushing), hence, we certainly pushed at least $\bm
			c_{\min} \leq \bm a_{i^*-1,w}$ into the first arc $v^*w$ of $P$. (This is also true if $v^* = s$.) For the remaining
			arcs of $P$ it is true, because by the lexicographical minimality of $P$, the algorithm always pushes a flow value that is greater or equal to $\bm c_{\min}$ first along the next arc on $P$.
			
			During the clean-up, this property remains valid as we only reduce flow on arcs that have a distance of more than $i^*$ from $t$.
			
			Hence, an arc $e \in P$ with $\bm c_e = \bm c_{\min}$ is saturated at the very end of the subroutine.
		\end{subproof}
		In conclusion, $\bm y$ is a feasible $s$-$t$-flow in the residual network that has positive value only on paths of length $k$ and saturates at least one arc. This finalizes the proof of \Cref{Thm:augmenting_flow}.
	\end{proof}

	\mypara{Correctness of the Main Routine}
	
	The following theorem states that \Cref{Alg:max_flow} correctly computes a maximum flow. Using the correctness of the technically involved subroutine \texttt{FindAugmentingFlow$_k$}, the remaining proof is actually similar to textbook proofs for the algorithms by Edmonds-Karp and Dinic; see, e.g., \cite{kortevygen}.
	
	\begin{restatable}{theorem}{thmmaapcorrectness}\label{Thm:MAAPCorrectness}
		Let $G=(V,E)$ be a fixed directed graph with $s,t\in V$. For capacities $(\bm\nu_e)_{e \in E}$ as input, the MAAP given by \Cref{Alg:max_flow} returns a maximum $s$-$t$-flow~$(\bm x_{e})_{e \in \vec E}$.
	\end{restatable}
	\begin{proof}
		It is a well-known fact that a feasible $s$-$t$-flow is maximum if and only if the corresponding residual network does not contain any $s$-$t$-path, see e.g.~\cite[Theorem~8.5]{kortevygen}. Since any simple path has length at most $n-1$, it suffices to show the following claim.
		\begin{claim}\label{Claim:outer}
			After iteration~$k$ of the \textbf{for} loop in line~\ref{Line:outerFor} of \Cref{Alg:max_flow}, $\bm x$ is a feasible $s$-$t$-flow with corresponding residual capacities $\bm c$ such that no $s$-$t$-path of length at most $k$ remains in the residual network.
		\end{claim}
		
		Given a residual network $(V,E^*)$, let $E^*_k$ be the set of arcs that lie on an $s$-$t$-path of length exactly $k$. If the distance from $s$ to $t$ is exactly $k$, then these arcs coincide with the arcs of the so-called \emph{level graph} used in Dinic's algorithm, compare \cite{dinic1970algorithm, kortevygen}.
		
		We will show \Cref{Claim:outer} about the outer \textbf{for} loop by induction on $k$ using a similar claim about the inner \textbf{for} loop.	
		
		\begin{claim}\label{Claim:inner}
			Suppose, at the beginning of an iteration of the \textbf{for} loop in line~\ref{Line:innerFor}, it holds that
			\begin{enumerate}[(i)]
				\item $\bm x$ is a feasible $s$-$t$-flow with corresponding residual capacities $\bm c$, and\label{firstitem}
				\item the length of the shortest $s$-$t$-path in the residual network is at least $k$.\label{lastitem}
			\end{enumerate}
			Then, after that iteration, properties \eqref{firstitem} and \eqref{lastitem} do still hold. Moreover, if $E^*_k$ is nonempty, then its cardinality is strictly reduced by that iteration.
		\end{claim}
		
		\begin{subproof}[Proof of \Cref{Claim:inner}]
			Since \eqref{firstitem} and $\eqref{lastitem}$ hold at the beginning of the iteration, \Cref{Thm:augmenting_flow} implies that the flow $\bm y$ found in line~\ref{Line:sub} fulfills flow conservation and is bounded by $-\bm c_{vu} \leq \bm y_{uv} \leq \bm c_{uv}$ for each $uv\in \vec{E}$. Hence, we obtain that, after updating $\bm x$ and $\bm c$ in lines \ref{Line:startUpdate} to \ref{Line:endUpdate}, $\bm x$ is still a feasible flow that respects flow conservation and capacities, and $\bm c$ are the corresponding new residual capacities. Thus, property \eqref{firstitem} is also true at the end of the iteration.
			
			Let $G^*=(V,E^*)$ and $\tilde{G}^*=(V,\tilde{E}^*)$ be the residual graphs before and after the iteration, respectively. Let $E^*_k$ and $\tilde E^*_k$ be the set of arcs on $s$-$t$-paths of length $k$ in $G^*$ and $\tilde{G}^*$, respectively. Finally, let $E'$ be the union of $E^*$ with the reverse arcs of $E^*_k$ and let $G'=(V,E')$.
			
			Since, by \Cref{Thm:augmenting_flow}, we only augment along arcs in $E^*_k$, it follows that $\tilde{E}^*\subseteq E'$. Let $P$ be a shortest $s$-$t$-path in $G'$ and suppose for contradiction that $P$ contains an arc that is not in $E^*$. Let $e=uv$ be the first of all such arcs of $P$ and let $P_u$ be the subpath of $P$ until node $u$. Then the reverse arc $vu$ must be in $E^*_k$. In particular, $\dist_{G^*}(s,v)<\dist_{G^*}(s,u)\leq \abs{E(P_u)}$, where the second inequality follows because $P_u$ uses only arcs in $E^*$. Hence, replacing the part of $P$ from $s$ to $v$ by a shortest $s$-$v$-path in $G^*$ reduces the length of $P$ by at least two, contradicting that $P$ is a shortest path in~$G'$.
			
			Thus, all shortest paths in $G'$ only contain arcs from $E^*$. In particular, they have length at least $k$. Hence, all paths in $G'$ that contain an arc that is not in $E^*$ have length larger than $k$. Since $\tilde{E}^*\subseteq E'$, this also holds for paths in $\tilde{G}^*$, which implies \eqref{lastitem}. It also implies that $\tilde E^*_k\subseteq E^*_k$. Moreover, by \Cref{Thm:augmenting_flow}, if $E^*_k$ is nonempty, at least one arc of $E^*_k$ is saturated during the iteration, and thus removed from $E^*_k$. Thus, the cardinality of $E^*_k$ becomes strictly smaller.
		\end{subproof}
		
		Using \Cref{Claim:inner}, we are now able to show \Cref{Claim:outer}.
		
		\begin{subproof}[Proof of \Cref{Claim:outer}]
			We use induction on~$k$. For the induction start, note that before entering the \textbf{for} loop in line~\ref{Line:outerFor}, that is, so to speak, after iteration $0$, obviously no $s$-$t$-path of length~$0$ can exist in the residual network. Also note that after the initialization in lines~\ref{Line:startIni} to~\ref{Line:endIni}, $\bm x$ is the zero flow, which is obviously feasible, and $\bm c$ contains the corresponding residual capacities.
			
			For the induction step, consider the $k$-th iteration. By the induction hypothesis, we know that, at the beginning of the $k$-th iteration, $\bm x$ is a feasible $s$-$t$-flow with corresponding residual capacities $\bm c$ and the distance from $s$ to $t$ in the residual network is at least $k$. Observe that by \Cref{Claim:inner}, these properties are maintained throughout the entire $k$-th iteration. In addition, observe that at the beginning of the $k$-th iteration, we have $\abs{E^*_k}\leq m$. Since, due to \Cref{Claim:inner}, $\abs{E^*_k}$ strictly decreases with each inner iteration until it is zero, it follows that after the $m$ inner iterations, the residual network does not contain an $s$-$t$-path of length $k$ any more, which completes the induction.
		\end{subproof}
		
		Since any simple path has length at most~$n-1$, \Cref{Claim:outer} implies that, at the end of iteration~$k=n-1$, the nodes $s$ and $t$ must be disconnected in the residual network. Hence, \Cref{Alg:max_flow} returns a maximum flow, which concludes the proof of \Cref{Thm:MAAPCorrectness}.
	\end{proof}

	\mypara{Bounding the Complexity}
	
	 By applying the definition of our complexity measures to \Cref{Alg:max_flow} and the subroutine, we first prove the following bounds.
	
	\begin{restatable}{theorem}{thmcomplexity}\label{Thm:MAAPComplexity}
		For the complexity measures of the MAAP $A$ defined by \Cref{Alg:max_flow} it holds that $d(A), s(A) \in\mathcal{O}(n^2m^2)$ and~\mbox{$w(A) \in \mathcal{O}(n^2)$}.
	\end{restatable}
	\begin{proof}
		We first analyze the MAAP~$A'$ given in \Cref{Alg:augmenting_flow}.	
		Concerning the complexity measures $d$ and $s$, the bottleneck of \Cref{Alg:augmenting_flow} is given by the two blocks consisting of lines~\ref{line:outer_for} to~\ref{line:increase_of_Y} as well as lines~\ref{line:begin_clean_up} to~\ref{line:increase_of_Y_by_b}. Each of these blocks has $\mathcal{O}(km)$ sequential iterations and the body of the innermost \textbf{for} loop has constant complexity. Thus, we have $d(A'), s(A')\in\mathcal{O}(km)\subseteq\mathcal{O}(nm)$ for the overall subroutine.
		
		Concerning measure $w$, the bottleneck is in fact the initialization in lines~\ref{line:start} to~\ref{line:init_fattest}, such that we have $w(A')\in\mathcal{O}(m+kn)\subseteq\mathcal{O}(n^2)$.
		
		Now consider the main routine in \Cref{Alg:max_flow}. For all three complexity measures, the bottleneck is the call of the subroutine in line~\ref{Line:sub} within the two \textbf{for} loops. Since we have a total of $\mathcal{O}(nm)$ sequential iterations, the claimed complexity measures follow. Note that the parallel \textbf{for} loops in lines~\ref{Line:startIni} and~\ref{Line:startUpdate} do not increase the measure $w(A)\in\mathcal{O}(n^2)$.
	\end{proof}
	
	We remark that the reported complexity $w(A)\in\mathcal{O}(n^2)$ in \Cref{Thm:MAAPComplexity} is actually suboptimal and can be replaced with $\mathcal{O}(1)$ instead, if the parallel \textbf{for} loops in the MAAP and the subroutine are replaced with sequential ones. The asymptotics of $d(A)$ and $s(A)$ remain unchanged because the bottleneck parts of the MAAP are already of sequential nature. Still, we used parallel \textbf{for} loops in \Cref{Alg:max_flow} and the subroutine in order to point out at which points the ability of NNs to parallelize can be used in a straightforward way.
	Combining the previous observations with \Cref{Prop:NN-MAAP} we obtain the following corollary.
	
	\begin{corollary}
		Let $G=(V,E)$ be a fixed directed graph with $s,t\in V$, $\abs{V}=n$, and $\abs{E}=m$. There exists an NN of depth and size $\mathcal{O}(m^2n^2)$ and width $\mathcal{O}(1)$ that correctly maps arc capacities $(\bm \nu_e)_{e \in E}$ to a maximum $s$-$t$-flow $(\bm x_{e})_{e \in \vec E}$.
	\end{corollary}
	
	Note that the objective value can be computed from the solution by simply adding up all flow values of arcs leaving the source node $s$. Therefore, this can be done by an NN with the same asymptotic size bounds, too.
	
	Finally, observe that the total computational work of \Cref{Alg:max_flow}, represented by $s(A)$ and, equivalently, the size of the resulting NN, differs only by a factor of $n$ from the standard running time bound $\mathcal{O}(nm^2)$ of the Edmonds-Karp algorithm; see~\cite[Corollary~8.15]{kortevygen}. While the number of augmenting steps is in $\mathcal{O}(nm)$ for both algorithms, the difference lies in finding the augmenting flow. While the Edmonds-Karp algorithm finds the shortest path in the residual network in $\mathcal{O}(m)$ time, our subroutine requires $\mathcal{O}(mn)$ computational work.
	
	\section{Future Research}
	
	Our grand vision, to which we contribute in this paper, is to determine how the computational model defined by ReLU~NNs compares to strongly polynomial time algorithms. In other words, is there a CPWL function (related to a CO problem or not) which can be evaluated in strongly polynomial time, but for which no polyonmial-size NNs exist?
	Resolving this question involves, of course, the probably very challenging task to prove nontrivial lower bounds on the size of NNs to compute certain functions. Particular candidate problems, for which the existence of polynomial-size NNs is open, are, for example, the Assignment Problem, different versions of Weighted Matching Problems, or Minimum Cost Flow Problems.
	
	Another direct question for further research is to what extend the size of the NN constructions in this paper can be improved. For example, is the size of $\mathcal{O}(n^2m^2)$ to compute maximum flows best possible or are smaller constructions conceivable? Even though highly parallel architectures (polylogarithmic depth) with polynomial width are unlikely,
	is it still possible to make use of NNs' ability to parallelize and find a construction with a depth that is a polynomial of lower degree than $n^2m^2$?
	
	We hope that this paper promotes further research on the interplay between NNs and classical (combinatorial optimization) algorithms.
	
	
	
	
	\bibliography{neuralmaxflow}


\end{document}

neuralmaxflow.bbl
\begin{thebibliography}{10}

\bibitem{AhujaMagnantiOrlin:NetworkFlows}
Ravindra~K. Ahuja, Thomas~L. Magnanti, and James~B. Orlin.
\newblock {\em Network flows: theory, algorithms, and applications}.
\newblock Prentice Hall, Upper Saddle River, New Jersey, USA, 1993.

\bibitem{ali1991neural}
M.K.~Mehmet Ali and Faouzi Kamoun.
\newblock A neural network approach to the maximum flow problem.
\newblock In {\em IEEE Global Telecommunications Conference GLOBECOM'91:
  Countdown to the New Millennium. Conference Record}, pages 130--134, 1991.

\bibitem{anthony2009neural}
Martin Anthony and Peter~L. Bartlett.
\newblock {\em Neural network learning: Theoretical foundations}.
\newblock Cambridge University Press, 1999.

\bibitem{Arora:DNNwithReLU}
Raman Arora, Amitabh Basu, Poorya Mianjy, and Anirbit Mukherjee.
\newblock Understanding deep neural networks with rectified linear units.
\newblock In {\em International Conference on Learning Representations}, 2018.

\bibitem{arora2009computational}
Sanjeev Arora and Boaz Barak.
\newblock {\em Computational complexity: a modern approach}.
\newblock Cambridge University Press, 2009.

\bibitem{beiu1996circuit}
Valeriu Beiu and John~G Taylor.
\newblock On the circuit complexity of sigmoid feedforward neural networks.
\newblock {\em Neural Networks}, 9(7):1155--1171, 1996.

\bibitem{Bello:NeuralCOwithRL}
Irwan {Bello}, Hieu {Pham}, Quoc~V. {Le}, Mohammad {Norouzi}, and Samy
  {Bengio}.
\newblock Neural combinatorial optimization with reinforcement learning.
\newblock {\em arXiv:1611.09940}, 2016.

\bibitem{BengioLodiProuvost:MLforCO}
Yoshua {Bengio}, Andrea {Lodi}, and Antoine {Prouvost}.
\newblock Machine learning for combinatorial optimization: a methodological
  tour d'horizon.
\newblock {\em arXiv:1811.06128}, 2018.

\bibitem{cybenko1989approximation}
George Cybenko.
\newblock Approximation by superpositions of a sigmoidal function.
\newblock {\em Mathematics of control, signals and systems}, 2(4):303--314,
  1989.

\bibitem{dinic1970algorithm}
Efim~A. Dinic.
\newblock Algorithm for solution of a problem of maximum flow in a network with
  power estimation.
\newblock {\em Soviet Mathematics Doklady}, 11:1277--1280, 1970.

\bibitem{edmonds1972theoretical}
Jack Edmonds and Richard~M. Karp.
\newblock Theoretical improvements in algorithmic efficiency for network flow
  problems.
\newblock {\em Journal of the ACM}, 19(2):248--264, 1972.

\bibitem{effati2008neural}
Sohrab Effati and M~Ranjbar.
\newblock Neural network models for solving the maximum flow problem.
\newblock {\em Applications and Applied Mathematics}, 3(3):149--162, 2008.

\bibitem{eldan2016power}
Ronen Eldan and Ohad Shamir.
\newblock The power of depth for feedforward neural networks.
\newblock In {\em Conference on Learning Theory}, pages 907--940, 2016.

\bibitem{Emami:learningPermutations}
Patrick {Emami} and Sanjay {Ranka}.
\newblock Learning permutations with sinkhorn policy gradient.
\newblock {\em arXiv:1805.07010}, 2018.

\bibitem{fomin2016subtraction}
Sergey Fomin, Dima Grigoriev, and Gleb Koshevoy.
\newblock Subtraction-free complexity, cluster transformations, and spanning
  trees.
\newblock {\em Foundations of Computational Mathematics}, 16(1):1--31, 2016.

\bibitem{gallo1989fast}
Giorgio Gallo, Michael~D Grigoriadis, and Robert~E Tarjan.
\newblock A fast parametric maximum flow algorithm and applications.
\newblock {\em SIAM Journal on Computing}, 18(1):30--55, 1989.

\bibitem{glorot2011deep}
Xavier Glorot, Antoine Bordes, and Yoshua Bengio.
\newblock Deep sparse rectifier neural networks.
\newblock In {\em Proceedings of the fourteenth international conference on
  artificial intelligence and statistics}, pages 315--323, 2011.

\bibitem{goldschlager1982}
Leslie~M. Goldschlager, Ralph~A. Shaw, and John Staples.
\newblock The maximum flow problem is log space complete for {P}.
\newblock {\em Theoretical Computer Science}, 21(1):105--111, 1982.

\bibitem{greenlaw1995limits}
Raymond Greenlaw, H.~James Hoover, and Walter~L. Ruzzo.
\newblock {\em Limits to parallel computation: {P}-completeness theory}.
\newblock Oxford University Press, 1995.

\bibitem{hanin2019universal}
Boris Hanin.
\newblock Universal function approximation by deep neural nets with bounded
  width and {ReLU} activations.
\newblock {\em Mathematics}, 7(10):992, 2019.

\bibitem{hanin2019complexity}
Boris Hanin and David Rolnick.
\newblock Complexity of linear regions in deep networks.
\newblock In {\em International Conference on Machine Learning}, pages
  2596--2604, 2019.

\bibitem{hanin2017approximating}
Boris Hanin and Mark Sellke.
\newblock Approximating continuous functions by {ReLU} nets of minimal width.
\newblock {\em arXiv:1710.11278}, 2017.

\bibitem{hertrich2021towards}
Christoph Hertrich, Amitabh Basu, Marco Di~Summa, and Martin Skutella.
\newblock Towards lower bounds on the depth of {ReLU} neural networks.
\newblock In {\em Proceedings of the 35th Conference on Neural Information
  Processing Systems (to appear)}, 2021.
\newblock Full version: \emph{arXiv:2105.14835}.

\bibitem{knapsackPaper}
Christoph Hertrich and Martin Skutella.
\newblock Provably good solutions to the knapsack problem via neural networks
  of bounded size.
\newblock In {\em Thirty-Fifth AAAI Conference on Artificial Intelligence (to
  appear)}, 2021.
\newblock Full version: \emph{arXiv:2005.14105}.

\bibitem{hopfield1985neural}
John~J. Hopfield and David~W. Tank.
\newblock ``{N}eural'' computation of decisions in optimization problems.
\newblock {\em Biological Cybernetics}, 52(3):141--152, 1985.

\bibitem{hornik1991approximation}
Kurt Hornik.
\newblock Approximation capabilities of multilayer feedforward networks.
\newblock {\em Neural networks}, 4(2):251--257, 1991.

\bibitem{jerrum1982some}
Mark Jerrum and Marc Snir.
\newblock Some exact complexity results for straight-line computations over
  semirings.
\newblock {\em Journal of the ACM (JACM)}, 29(3):874--897, 1982.

\bibitem{jukna2015lower}
Stasys Jukna.
\newblock Lower bounds for tropical circuits and dynamic programs.
\newblock {\em Theory of Computing Systems}, 57(1):160--194, 2015.

\bibitem{jukna2019greedy}
Stasys Jukna and Hannes Seiwert.
\newblock Greedy can beat pure dynamic programming.
\newblock {\em Information Processing Letters}, 142:90--95, 2019.

\bibitem{kennedy1988neural}
Michael~Peter Kennedy and Leon~O Chua.
\newblock Neural networks for nonlinear programming.
\newblock {\em IEEE Transactions on Circuits and Systems}, 35(5):554--562,
  1988.

\bibitem{Khalil:LearningCOoverGraphs}
Elias Khalil, Hanjun Dai, Yuyu Zhang, Bistra Dilkina, and Le~Song.
\newblock Learning combinatorial optimization algorithms over graphs.
\newblock In I.~Guyon, U.~V. Luxburg, S.~Bengio, H.~Wallach, R.~Fergus,
  S.~Vishwanathan, and R.~Garnett, editors, {\em Advances in Neural Information
  Processing Systems 30}, pages 6348--6358. 2017.

\bibitem{kool2019attention}
Wouter Kool, Herke van Hoof, and Max Welling.
\newblock Attention, learn to solve routing problems!
\newblock In {\em International Conference on Learning Representations}, 2019.

\bibitem{kortevygen}
Bernhard Korte and Jens Vygen.
\newblock {\em Combinatorial Optimization: Theory and Algorithms}.
\newblock Springer, 4th edition, 2008.

\bibitem{LeCunBengioHinton:DeepLearning}
Yann {LeCun}, Yoshua Bengio, and Geoffrey Hinton.
\newblock Deep learning.
\newblock {\em Nature}, 521:436--444, 2015.

\bibitem{liang2017deep}
Shiyu Liang and Rayadurgam Srikant.
\newblock Why deep neural networks for function approximation?
\newblock In {\em International Conference on Learning Representations}, 2017.

\bibitem{lodi2017learning}
Andrea Lodi and Giulia Zarpellon.
\newblock On learning and branching: a survey.
\newblock {\em {TOP}}, 25(2):207--236, 2017.

\bibitem{mccormick1999fast}
S~Thomas McCormick.
\newblock Fast algorithms for parametric scheduling come from extensions to
  parametric maximum flow.
\newblock {\em Operations Research}, 47(5):744--756, 1999.

\bibitem{montufar2014regions}
Guido~F. Montufar, Razvan Pascanu, Kyunghyun Cho, and Yoshua Bengio.
\newblock On the number of linear regions of deep neural networks.
\newblock In Z.~Ghahramani, M.~Welling, C.~Cortes, N.~D. Lawrence, and K.~Q.
  Weinberger, editors, {\em Advances in Neural Information Processing Systems
  27}, pages 2924--2932. 2014.

\bibitem{mukherjee2017lower}
Anirbit Mukherjee and Amitabh Basu.
\newblock Lower bounds over boolean inputs for deep neural networks with {ReLU}
  gates.
\newblock {\em arXiv:1711.03073}, 2017.

\bibitem{nazemi2012}
Alireza Nazemi and Farahnaz Omidi.
\newblock A capable neural network model for solving the maximum flow problem.
\newblock {\em Journal of Computational and Applied Mathematics},
  236(14):3498--3513, 2012.

\bibitem{nguyen2018neural}
Quynh Nguyen, Mahesh~Chandra Mukkamala, and Matthias Hein.
\newblock Neural networks should be wide enough to learn disconnected decision
  regions.
\newblock In {\em International Conference on Machine Learning}, pages
  3737--3746, 2018.

\bibitem{nowak:quadrAssignment}
Alex {Nowak}, Soledad {Villar}, Afonso~S. {Bandeira}, and Joan {Bruna}.
\newblock {Revised Note on Learning Algorithms for Quadratic Assignment with
  Graph Neural Networks}.
\newblock {\em arXiv:1706.07450}, 2017.

\bibitem{orlin2013}
James~B. Orlin.
\newblock Max flows in {O}(nm) time, or better.
\newblock In {\em Proceedings of the Forty-Fifth Annual ACM Symposium on Theory
  of Computing (STOC '13)}, pages 765--774. Association for Computing
  Machinery, 2013.

\bibitem{parberry1994circuit}
Ian Parberry, Michael~R Garey, and Albert Meyer.
\newblock {\em Circuit complexity and neural networks}.
\newblock MIT Press, 1994.

\bibitem{pascanu2014number}
Razvan Pascanu, Guido Montufar, and Yoshua Bengio.
\newblock On the number of inference regions of deep feed forward networks with
  piece-wise linear activations.
\newblock In {\em International Conference on Learning Representations}, 2014.

\bibitem{raghu2017expressive}
Maithra Raghu, Ben Poole, Jon Kleinberg, Surya Ganguli, and Jascha~Sohl
  Dickstein.
\newblock On the expressive power of deep neural networks.
\newblock In {\em International Conference on Machine Learning}, pages
  2847--2854, 2017.

\bibitem{safran2017depth}
Itay Safran and Ohad Shamir.
\newblock Depth-width tradeoffs in approximating natural functions with neural
  networks.
\newblock In {\em International Conference on Machine Learning}, pages
  2979--2987, 2017.

\bibitem{serra2018bounding}
Thiago Serra, Christian Tjandraatmadja, and Srikumar Ramalingam.
\newblock Bounding and counting linear regions of deep neural networks.
\newblock In {\em International Conference on Machine Learning}, pages
  4565--4573, 2018.

\bibitem{Shalev2014:UnderstandingML}
Shai Shalev-Shwartz and Shai Ben-David.
\newblock {\em Understanding machine learning: From theory to algorithms}.
\newblock Cambridge University Press, 2014.

\bibitem{shawe1992classes}
John~S Shawe-Taylor, Martin~HG Anthony, and Walter Kern.
\newblock Classes of feedforward neural networks and their circuit complexity.
\newblock {\em Neural networks}, 5(6):971--977, 1992.

\bibitem{shpilka2010arithmetic}
Amir Shpilka and Amir Yehudayoff.
\newblock {\em Arithmetic circuits: A survey of recent results and open
  questions}.
\newblock Now Publishers Inc, 2010.

\bibitem{Smith:NNforCOreview}
Kate~A. Smith.
\newblock Neural networks for combinatorial optimization: A review of more than
  a decade of research.
\newblock {\em INFORMS Journal on Computing}, 11(1):15--34, 1999.

\bibitem{Telgarsky15}
Matus Telgarsky.
\newblock Representation benefits of deep feedforward networks.
\newblock {\em arXiv:1509.08101}, 2015.

\bibitem{telgarsky2016benefits}
Matus Telgarsky.
\newblock Benefits of depth in neural networks.
\newblock In {\em Conference on Learning Theory}, pages 1517--1539, 2016.

\bibitem{Vinyals:PointerNetworks}
Oriol Vinyals, Meire Fortunato, and Navdeep Jaitly.
\newblock Pointer networks.
\newblock In C.~Cortes, N.~D. Lawrence, D.~D. Lee, M.~Sugiyama, and R.~Garnett,
  editors, {\em Advances in Neural Information Processing Systems 28}, pages
  2692--2700. 2015.

\bibitem{williamson_2019}
David~P. Williamson.
\newblock {\em Network Flow Algorithms}.
\newblock Cambridge University Press, 2019.

\bibitem{yarotsky2017error}
Dmitry Yarotsky.
\newblock Error bounds for approximations with deep relu networks.
\newblock {\em Neural Networks}, 94:103--114, 2017.

\end{thebibliography}
