# Breaking the $n$-Pass Barrier: A Streaming Algorithm for Maximum Weight
Bipartite Matching

S. Cliff Liu Princeton University<EMAIL_ADDRESS>Zhao Song School of
Mathematics, Institute for Advanced Study<EMAIL_ADDRESS>Hengjie Zhang
Columbia University<EMAIL_ADDRESS>
Given a weighted bipartite graph with $n$ vertices and $m$ edges, the _maximum
weight bipartite matching_ problem is to find a set of vertex-disjoint edges
with the maximum weight. This classic problem has been extensively studied for
over a century.

In this paper, we present a new streaming algorithm for the maximum weight
bipartite matching problem that uses $\widetilde{O}(n)$ space and
$\widetilde{O}(\sqrt{m})$ passes, which breaks the $n$-pass barrier. All the
previous streaming algorithms either require $\Omega(n\log n)$ passes or only
find an approximate solution. Our streaming algorithm constructs a subgraph
with $n$ edges of the input graph in $\widetilde{O}(\sqrt{m})$ passes, such
that the subgraph admits the optimal matching with good probability.

Our method combines various ideas from different fields, most notably the
construction of _space-efficient_ interior point method (IPM), SDD system
solvers, the isolation lemma, and LP duality. To the best of our knowledge,
this is the first work that implements the SDD solvers and IPMs in the
streaming model in $\widetilde{O}(n)$ spaces for graph matrices; previous IPM
algorithms only focus on optimizing the running time, regardless of the space
usage.

##  1 Introduction

The growing need for processing massive graphs in modern applications has
crowned the _data streaming_ model as a desirable testbed for scalable
algorithm design for graph problems. Since the introduction of the graph
streaming model in [FKM+04], there has been a surging line of work studying
basic graph problems in this model, such as connectivity [AGM12a], max-cut
[KKS14b, KK15, KKSV17], spectral sparsifiers [KL11, AGM12b, KLM+17, KMM+19,
KNST19], spanner construction [AGM12b, KW14], and so on.

Bipartite matching, which is the focus of this work, is one of the most
fundamental problems in theoretical computer science, dating back over a
century. Broadly speaking, there have been two conceptual approaches to
designing efficient matching algorithms – the traditional one is to leverage
the _combinatorial_ properties in the graph [Tar83]; the other is
reformulating the graph problem as a specific linear program [DS08], and using
_continuous optimization_ techniques (e.g., interior point method [Kar84] and
SDD system solver [ST04]) to design efficient algorithms.

In the streaming setting, the _bipartite matching_ problem has been
extensively studied [FKM+04, AG11, EKMS12, GKK12, EKMS12, Kap13, GO13, DNO14,
KKS14a, DNO14, AG18, ABB+19, AKSY20, AR20, FHM+20, Ber20, CKP+21]. Many of
these works consider the _semi-streaming_ model, where each edge is revealed
along with its weight one-by-one in an online fashion and according to an
adversarial order, and the algorithm is allowed to make one or more passes
over the stream using $\widetilde{O}(n)$ space.111Some authors define the
space in the streaming model to be the number of cells, where each cell can
hold $O(\log n)$ bits or even a number with infinite precision. Our bounds
remain unchanged even if each cell only holds a constant number of bits. We
remark that the order of edges in different passes _need not_ to be the same.
Despite the intimate relationship between this model and the matching problem,
almost all existing streaming algorithms for maximum weight bipartite matching
are only able to find an _approximate_ solution. To the best of our knowledge,
the only exception is a folklore algorithm that repeatedly finds an augmenting
path by a breath-first search, which takes $O(n\log n)$ passes and
$\widetilde{O}(n)$ space to find an exact maximum cardinality bipartite
matching. Obtaining even an $O(n)$-pass semi-streaming algorithm for this
problem seems to be a major challenge, which raises the following open
problem:

> _Is there a semi-streaming algorithm for exact maximum weight bipartite
> matching that uses $o(n)$ passes?_

Another motivation for this problem is the two (folklore) reductions: (1)
$s$-$t$ reachability $\preceq$ maximum bipartite matching, and (2) $s$-$t$
reachability $\preceq$ single source shortest path, all of which are
fundamental problems in algorithmic graph theory. We note that $s$-$t$
reachability [LJS19] and single source shortest path [FKM+09, KM19, CFCHT20]
can be both solved in $n^{1/2+o(1)}$ passes in semi-streaming, whereas for
exact maximum cardinality bipartite matching, only an $O(n\log n)$ passes
upper bound is known, leaving a huge open gap. To the best of our knowledge,
[KM19, CFCHT20] are also the only two semi-streaming algorithms with an
$n^{\Omega(1)}$ pass bounds. As such, a further motivation for answering our
problem is a better understanding on the power of semi-streaming in the multi-
pass regime.

Our main result is an affirmative answer to the above problem, obtaining a
sublinear-pass bipartite matching algorithm through _space-efficient interior
point method_. The interior point method (IPM) optimization framework was
proposed by Karmarkar [Kar84] as a powerful tool for improving the running
time of linear programming solver and, subsequently, for several specific
graph-type linear programming. Alas, in contrast to classical IPM solvers,
which do not have a space constraint, in the streaming setting the algorithm
is not allowed to store all the $m$ constraints of the ($n$-variable) LP.

Indeed, the most enduring message of our result is showing that _interior
point method for maximum bipartite matching can be implemented in the semi-
streaming model in a non-trivial space efficient manner, whereas for general
matrices, LP solvers are impossible to be implemented in near-linear space_.
More precisely, we present a semi-streaming algorithm that takes
$\widetilde{O}(\sqrt{m})$ IPM iterations to compute the maximum weight
bipartite matching, where each iteration employs $O(\log n)$ passes to
implement a streaming version of an SDD solver for computing the Newton step
of the path-following IPM. To the best of our knowledge, our algorithm is the
first to implement IPM in space-limited models such as the streaming model.

More specifically, to solve the maximum weight bipartite matching problem in
$\widetilde{O}(n)$ space, we consider its dual, the _generalized minimum
vertex cover_ problem in bipartite graph, and solve this problem using a
streaming version of IPM. Therefore, this work also yields an
$\widetilde{O}(\sqrt{m})$-pass semi-streaming algorithm for the minimum vertex
cover problem in bipartite graph. To transform the dual solution to a primal
solution in $\widetilde{O}(n)$ space, we use the isolation lemma, whose
original construction uses $\widetilde{\Theta}(m)$ random bits [MVV87]. We
bypass this issue by implementing an alternative construction proposed by
[CRS95] that uses $\widetilde{O}(n)$ random bits in the semi-streaming model.

###  1.1 Our contribution

We give our main result and two byproducts as follows.

We provide a semi-streaming algorithm for bipartite matching running in
$\widetilde{O}(\sqrt{m})$ passes. This solves the long standing open problem
of whether maximum matching can be solved in $o(n)$ passes, in all regimes
where $m=n^{2-\epsilon}$ for any constant $\epsilon>0$.

######  Theorem 1.1 (Main theorem, informal version of Theorem 11.1).

Given a bipartite graph $G$ with $n$ vertices and $m$ edges, there exists a
streaming algorithm that computes an (exact) maximum weight matching of $G$ in
$\widetilde{O}(\sqrt{m})$ passes and $\widetilde{O}(n)$ space with probability
$1-1/\operatorname{poly}(n)$.

We also show that minimum vertex cover in bipartite graph, which is the dual
problem of maximum matching, can be solved within the same number of passes
and space. A _vertex cover_ of a graph is a set of vertices that includes at
least one endpoint of every edge of the graph, and a _minimum vertex cover_ is
a vertex cover of minimum size.

######  Theorem 1.2 (Informal version of Theorem 10.8).

Given a bipartite graph $G$ with $n$ vertices and $m$ edges, there exists a
streaming algorithm that computes a minimum vertex cover of $G$ in
$\widetilde{O}(\sqrt{m})$ passes and $\widetilde{O}(n)$ space with probability
$1-1/\operatorname{poly}(n)$.222We can actually solve a _generalized_ version
of the minimum vertex cover problem in bipartite graph: each edge $e$ needs to
be _covered_ for at least $b_{e}\in\operatorname*{\mathbb{Z}}^{+}$ times,
where the case of $b=\mathbf{1}_{m}$ is the classic minimum vertex cover.

Our main result is obtained by a novel application of IPM together with SDD
system solver in the streaming model. To the best of our knowledge, this is
the first work that implements those techniques in the streaming model in
which only $\widetilde{O}(n)$ space is allowed. We put our result of SDD
solver here for possible independent interests.

######  Theorem 1.3 (Informal version of Theorem 9.9).

There is a streaming algorithm which takes input an SDD matrix $A$, a vector
$b\in\mathbb{R}^{n}$, and a parameter $\epsilon\in(0,1)$, if there exists
$x^{*}\in\mathbb{R}^{n}$ such that $Ax^{*}=b$, then with probability
$1-1/\operatorname{poly}(n)$, the algorithm returns an $x\in\mathbb{R}^{n}$
such that

$\|x-x^{*}\|_{A}\leq\epsilon\cdot\|x^{*}\|_{A}$

in $O(1+\log(1/\epsilon)/\log\log n)$ passes and $\widetilde{O}(n)$
space.333Our algorithm can actually solve any SDD0 system, which is more
general than an SDD system. See a formal definition of SDD0 matrix in Section
3.

We hope the tools developed in this paper will find further applications in
the study of streaming algorithms and beyond.

###  1.2 Related work

##### Streaming algorithms for matching.

Maximum matching has been extensively studied in the streaming model for
decades, where almost all of them fall into the category of approximation
algorithms. For algorithms that only make one pass over the edges stream,
researchers make continuous progress on pushing the constant approximation
ratio above $1/2$, which is under the assumption that the edges are arrived in
a uniform random order [KKS14a, ABB+19, FHM+20, Ber20]. The random-order
assumption makes the problem easier (at least algorithmically). A more general
setting is multi-pass streaming with adversarial edge arriving. Under this
setting, the first streaming algorithm that beats the $1/2$-approximation of
bipartite cardinality matching is [FKM+04], giving a
$2/3\cdot(1-\epsilon)$-approximation in $1/\epsilon\cdot\log(1/\epsilon)$
passes. The first to achieve a $(1-\epsilon)$-approximation is [McG05], which
takes $(1/\epsilon)^{1/\epsilon}$ passes.444For the weighted case, there is a
$(1/2-\epsilon)$-approximation algorithm that only takes one pass [PS17].
Since then, there is a long line of research in proving upper bounds and lower
bounds on the number of passes to compute a maximum matching in the streaming
model [AG11, EKMS12, GKK12, EKMS12, Kap13, DNO14, AG18, AKSY20, AR20, ALT21]
(see next subsection for more details). Notably, [AG11, AG18] use linear
programming and duality theory (see the next subsection for more details).

However, all the algorithms above can only compute an approximate maximum
matching: to compute a matching whose size is at least $(1-\epsilon)$ times
the optimal, one needs to spend $\text{poly}(1/\epsilon)$ passes (see [DNO14,
AG18] and the references therein).

Despite the intimate relationship between the matching problem and the
streaming model, no breakthrough is made towards solving exact bipartite
matching in the streaming model. We remark that a simple folklore algorithm
inspired by the classic algorithm of Hopcroft and Karp [HK73] can actually
find the exact maximum cardinality bipartite matching in $\widetilde{O}(n)$
passes using $\widetilde{O}(n)$ space, whose details can be found in the next
subsection.

##### Streaming spectral sparsifer.

Initialized by the study of cut sparsifier in the streaming model [AG09], a
simple one-pass semi-streaming algorithm for computing a spectral sparsifier
of any weighted graph is given in [KL11], which suffices for our use in this
paper. The problem becomes more challenging in a dynamic setting, i.e., both
insertion and deletion of edges from the graph are allowed. Using the idea of
linear sketching, [KLM+17] gives a single-pass semi-streaming algorithm for
computing the spectral sparsifier in the dynamic setting. However, their
brute-force approach to recover the sparsifier from the sketching uses
$\Omega(n^{2})$ time. An improved recover time is given in [KMM+19] but
requires more spaces, e.g., $\epsilon^{-2}n^{1.5}\log^{O(1)}n$. Finally,
[KNST19] proposes a single-pass semi-streaming algorithm that uses
$\epsilon^{-2}n\log^{O(1)}n$ space and $\epsilon^{-2}n\log^{O(1)}n$ recover
time to compute an $\epsilon$-spectral sparsifier which has
$O(\epsilon^{-2}n\log n)$ edges. Note that $\Omega(\epsilon^{-2}n\log n)$
space is necessary for this problem [CKST19].

##### SDD solver and IPM.

There is a long line of work focusing on fast SDD solvers [ST04, KMP10, KMP11,
KOSZ13, CKM+14, PS14, KS16]. Spielman and Teng give the first nearly-linear
time SDD solver, which is simplified with a better running time in later
works. The current fastest SDD solver runs in
$O(m\log^{1/2}n\operatorname{poly}(\log\log n)\log(1/\epsilon))$ time
[CKM+14]. All of them require $\widetilde{\Theta}(m)$ space.

The interior point method was originally proposed by Karmarkar [Kar84] for
solving linear program. Since then, there is a long line of work on speeding
up interior point method for solving classical optimization problems, e.g.,
linear program [Vai87, Ren88, Vai89, NN89, LS14, LS13, LS14, LS15, CLS19,
LSZ19, BLSS20, Bra20, SY21, JSWZ21, DLY21],555We remark that the fastest
classical LP solvers are due to [LS14], [JSWZ21], and [BLSS20]. The space
requirements for those algorithms are $O(nm)$ [LS14, BLSS20] and
$O((m+n)^{2})$ [JSWZ21]. and semi-definite program [Ans00, NN92, NN94, JKL+20,
HJST21].

###  1.3 Previous techniques

In this section, we summarize the previous techniques for computing the
maximum matching in the streaming model.

  * •

In Section 1.3.1, we introduce some representative approximation algorithms
for bipartite matching.

  * •

In Section 1.3.2, we present a potential method to compute an exact bipartite
matching, showcasing the current bottleneck in the field.

  * •

In Section 1.3.3, we discuss a simple folklore semi-streaming algorithm that
uses $O(n\log n)$ passes, which is by far the best.

  * •

In Section 1.3.4, we introduce some related bipartite matching algorithms
which use more space than that permitted in the semi-streaming model.

####  1.3.1 Approximation algorithms

Given a parameter $\epsilon\in(0,1)$, many streaming algorithms are to find a
matching of size $(1-\epsilon)$ times the size of the maximum matching. The
space and passes usages of these approximation algorithms are increasing
functions of $1/\epsilon$.666We will be focusing on approximate algorithms
that find a matching that is _close to_ (or can potentially be used to find)
an exact maximum matching, so all the constant-approximate algorithms are not
introduced here. We refer the interested readers to [AB19] and the references
therein.

A natural idea to find an approximate matching is to iteratively sample a
small subset of edges and use these edges to refine the current matching.
These algorithms are called _sampling-based_ algorithms. In [AG18], Ahn and
Guha show that by adaptively sampling $\widetilde{O}(n)$ edges in each
iteration, one can either obtain a certificate that the sampled edges admit a
desirable matching, or these edges can be used to refine the solution of a
specific LP. The LP is a nonstandard relaxation of the matching problem, and
will eventually be used to produce a good approximate matching. The algorithm
of Ahn and Guha can compute a $(1-\epsilon)$-approximate matching for weighted
(not necessarily bipartite) graph in $\widetilde{O}(1/\epsilon)$ passes and
$\widetilde{O}(n\operatorname{poly}(1/\epsilon))$ space. However, the degree
of $\operatorname{poly}(1/\epsilon)$ in the space usage can be very large,
making their algorithm inapplicable for small (non-constant)
$\epsilon=o(1/\log n)$ in the semi-streaming model.

Finding a $(1-\epsilon)$-approximate maximum matching with no space dependence
on $\epsilon$ requires different methods. Inspired by the well-studied _water
filling_ process in online algorithms (see [DJK13] and the references
therein), Kapralov proposes an algorithm that generalizes the water filling
process to multiple passes [Kap13]. This algorithm works in the vertex arrival
semi-streaming model, where a vertex and all of its incident edges arrive in
the stream together. The observation is that the water filling from pass
$(k-1)$ to pass $k$ follows the same manner as that in the first pass (with a
more careful double-counting method), then solving differential equations
gives a $(1-1/\sqrt{2\pi k})$-approximate matching in $k$ passes.

Kapralov’s algorithm removes the $\operatorname{poly}(\log n)$ factor in the
number of passes comparing to [AG11], giving a $(1-\epsilon)$-approximate
maximum matching in $O(1/\epsilon^{2})$ passes, albeit in a stronger vertex
arrival model. Very recently, Assadi, Liu, and Tarjan give a simple semi-
streaming algorithm that computes a $(1-\epsilon)$-approximate maximum
matching in $O(1/\epsilon^{2})$ passes, removing the vertex arrival condition
[ALT21]. Their algorithm considers the left vertices in a bipartite graph as
_bidders_ and the right vertices as _items_ , then the maximum cardinality
bipartite matching problem is to find an auction maximizing the social
welfare. The algorithm only needs to maintain a price for each item. In each
round, every unmatched bidder chooses the item that maximizes his own revenue.
They show that a simple auction scheme converges in $O(1/\epsilon^{2})$ rounds
(passes). The space usage is $O(n\log(1/\epsilon))$ for storing all the
prices, which is $\widetilde{O}(n)$ since $\epsilon\geq 1/(n+1)$ (otherwise we
could obtain an exact matching).

####  1.3.2 From approximate to exact maximum matching

A potential method to compute an exact maximum matching is to augment an
approximate matching by repeatedly finding augmenting paths.777Given a
matching in a graph, an _augmenting path_ is a path that starts and ends at an
unmatched vertex, and alternately contains edges that are outside and inside
the matching. The state-of-the-art semi-streaming algorithm based on auctions
by [ALT21] takes $O(1/\epsilon^{2})$ passes to find a
$(1-\epsilon)$-approximate maximum cardinality bipartite matching, where
$\epsilon\geq 1/(n+1)$ is a parameter. (For $\epsilon=\Omega(1/\log n)$, there
could be algorithms that take fewer passes. Recall that the algorithm of
[AG18] does not work in semi-streaming when $\epsilon$ is too small.) So given
$o(n)$ passes, no existing algorithm can find a matching of size at least
$\mathrm{OPT}-O(\sqrt{n})$, assuming the maximum matching has size
$\mathrm{OPT}=\Theta(n)$. Therefore, we are only able to obtain a matching of
size $\mathrm{OPT}-\Omega(\sqrt{n})$ from previous approximation algorithms.
However, currently there is no semi-streaming algorithm that solves directed
graph reachability – a problem that is no harder than finding one augmenting
path – in $o(\sqrt{n})$ passes [LJS19]. (The linear-work parallel algorithm of
[LJS19] can be translated into a semi-streaming algorithm that finds an
augmenting path in $n^{1/2+o(1)}$ passes.) So augmenting a matching of size
$\mathrm{OPT}-\Omega(\sqrt{n})$ to the maximum matching also takes $\Omega(n)$
passes. In conclusion, based on existing algorithms, the method of augmenting
an approximate matching does not yield an $o(n)$-pass semi-streaming
algorithm. To break the $n$-pass barrier under this framework, we need either
of the following:

  * •

An $o(n)$-pass approximation algorithm that finds a matching of size
$\mathrm{OPT}-n^{1/2-\Omega(1)}$.

  * •

An $n^{1/2-\Omega(1)}$-pass algorithm that finds at least one augmenting path.

####  1.3.3 A folklore algorithm with $O(n\log n)$ passes

A simple folklore algorithm inspired by the classic algorithm of Hopcroft and
Karp [HK73] can actually find the exact maximum cardinality bipartite matching
in $\widetilde{O}(n)$ passes using $\widetilde{O}(n)$ space. The main idea is
the following. Let $\mathrm{OPT}$ be the size of the maximum matching in the
given $n$-vertex bipartite graph. If the current matching has size $i$, then
there must exist $(\mathrm{OPT}-i)$ disjoint augmenting paths, so the shortest
augmenting path has length at most $n/(\mathrm{OPT}-i)$. Using a breath-first
search (simply ignore the edge in the stream that is not incident with the
frontier of the breath-first search), one can find this path in
$n/(\mathrm{OPT}-i)$ passes and augment the current matching. Therefore, the
total number of passes to compute the perfect matching is at most

$\sum_{i=0}^{\mathrm{OPT}-1}\frac{n}{\mathrm{OPT}-i}=n\cdot\sum_{i=1}^{\mathrm{OPT}}\frac{1}{i}=O(n\log
n).$

This simple algorithm is by far the state-of-the-art: to the best of our
knowledge, before this work, no (semi-)streaming algorithm computes an exact
matching in $o(n\log n)$ passes.

####  1.3.4 Other algorithms with more space

Since the breakthrough result from [Mad13], new improvements have been made in
solving the maximum weight bipartite matching problem, and more generally, the
maximum flow problem (see [LS20a] and the references therein). The maximum
cardinality bipartite matching algorithm by [LS20a] runs in $m^{4/3+o(1)}$
time and takes $\Omega(m)$ space, which is currently the fastest when the
graph is sufficiently sparse. For moderately dense graph, a better algorithm
is given in [BLN+20] very recently, which runs in $\widetilde{O}(m+n^{1.5})$
time and takes $\Omega(m)$ space. For the maximum weight bipartite matching
problem, the $O(n^{\omega})$-time algorithm using matrix inversion where
$\omega<2.373$ (see [MS04]) requires $\Omega(n^{2})$ space, and the
$\widetilde{O}(m\sqrt{n})$-time algorithm by [LS14] requires $\Omega(m)$
space.

We remark that, except IPM, another method to solve the matching problem is to
use simplex-type LP algorithms. However, those algorithms are mainly designed
for the case where the dimension $n$ is small, thus they do not produce
meaningful results under the graph setting where $n\leq m\leq n^{2}$. For
instance, [CC07] shows that the non-recursive implementation of Clarkson’s
algorithm [Cla95] gives a streaming LP solver that uses $O(n)$ passes and
$\widetilde{O}(n\sqrt{m})$ spaces, which is $\Omega(m)$ space for graph
problems. They also show that the recursive implementation gives a streaming
LP solver that uses $n^{O(1/\delta)}$ passes and
$(n^{2}+m^{\delta})\operatorname{poly}(1/\delta)$ spaces. Additionally,
[AKZ19] proposes a streaming algorithm for solving $n$-dimensional LP that
uses $O(nr)$ pass and $O(m^{1/r})\operatorname{poly}(n,\log m)$ space, where
$r\geq 1$ is a parameter.

##  2 Our techniques

Our approach is based on a novel combination of several techniques from both
continuous optimization and combinatorics, specifically, interior point
methods, streaming SDD solvers, LP duality, and the isolation lemma.

The key obstacle to overcome in this paper is the space limitation of
implementing IPMs in the streaming model.

To better appreciate the implication of the $\widetilde{O}(n)$ space
constraint, note that storing a matching already takes $\widetilde{\Theta}(n)$
space, meaning that we have only a polylog space overhead per vertex to store
auxiliary information. The conventional way of solving maximum bipartite
matching using an IPM solver would get stuck at the very beginning -
maintaining the solution of the relaxed linear program, which is a fractional
matching, already requires $\Omega(m)$ space for storing all LP constraints,
which seems inevitable.

Our main contribution is fixing this problem by showing that a certain IPM
([Ren88]) can actually be made to work in sublinear space (in the number of LP
constraints) for _graph_ constraint matrices. Our key insight is showing that
solving the _dual_ form of the above LP, which corresponds to the generalized
(fractional) minimum vertex cover problem, is sufficient, and therefore only
$\widetilde{O}(n)$ space is needed for maintaining a fractional solution. We
use several techniques to establish this argument. The first idea is to use
complementary slackness for the dual solution to learn which $n$ edges will be
in the final maximum matching and therefore reduce the size of the graph from
$m$ to $n$. However, this is not always the case: For instance, in a bipartite
graph that admits a perfect matching, all left vertices form a minimum vertex
cover, but the complementary slackness theorem gives no information on which
edges are in the perfect matching. To circumvent this, we need to slightly
perturb the weight on every edge, so that the minimum vertex cover (which is
now unique) indeed provides enough information. We use the isolation lemma
[MVV87] to achieve this.

Our second contribution is implementing the isolation lemma in limited space.
Perturbing the weight on every edge requires storing $\widetilde{O}(m)$ bits
of randomness, since the perturbation should remain identical across two
different passes. We bypass this issue by using the generalized isolation
lemma proposed by [CRS95], in which only $O(\log(Z))$ bits of randomness is
needed, where $Z$ is the number of candidates. In our case, $Z\leq n^{n}$ is
the number of all possible matchings. So $\widetilde{O}(n)$ space usage
perfectly fits into the semi-streaming model. We design an oracle that by
storing $\widetilde{O}(n)$ random bits, will output the same perturbations for
all edges in all passes.

Our third contribution is implementing an SDD (and Laplacian) solver in the
semi-streaming model, which, to the best of our knowledge, was not done prior
to our work.

The last piece of this work is to select a suitable IPM solver for our
framework. Many recently-developed fast IPM solvers do not fit into
$\widetilde{O}(n)$ space. For example, [LS20a, LS20b] cannot handle the dual
LP of matching, while the general LP solvers [LS14, JSWZ21, BLSS20] need to
maintain _both_ primal and dual solutions. The matching solver [BLN+20] based
on IPM also maintains both primal and dual solutions. Fortunately, we show
that Renegar’s IPM algorithm [Ren88] can be implemented into a streaming
version using all above techniques.

In the following subsections we elaborate on each of the above components:

  * •

In Section 2.1, we provide a high-level picture of how interior point method
works.

  * •

In Section 2.2, we describe our contribution on our implementations of SDD
solver, IPM, and the isolation lemma in the streaming model.

  * •

In Section 2.3, we show a novel application of the isolation lemma to turn
dual into primal.

###  2.1 IPM analysis

Given a matrix $A\in\mathbb{R}^{m\times n}$ and two vectors
$b\in\mathbb{R}^{m}$ and $c\in\mathbb{R}^{n}$, consider the following linear
program

$\displaystyle\min_{x\in\mathbb{R}^{n}}c^{\top}x\mathrm{\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ s.t.\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ }Ax\geq b.$

We briefly review several basic concepts in Renegar’s algorithm [Ren88] (see
Sections 5 and 6). Renegar’s IPM solves the above LP by introducing two
functions: one is the barrier function $F:\mathbb{R}^{n}\rightarrow\mathbb{R}$
and the other is the (parametrized) perturbed objective function
$f_{t}:\mathbb{R}^{n}\rightarrow\mathbb{R}$. The standard log-barrier function
can be written as

$F(x):=\sum_{i\in[m]}\log(a_{i}^{\top}x-b_{i}),$

where $a_{1}^{\top},\cdots,a_{m}^{\top}$ are rows of matrix $A$. It is easy to
see that the domain of $F(x)$ is exactly the region
$\\{x\in\mathbb{R}^{n}:Ax>b\\}$. For any parameter $t>0$, we define the
perturbed function $f_{t}:\mathbb{R}^{n}\rightarrow\mathbb{R}$ as follows

$\displaystyle f_{t}(x):=t\cdot c^{\top}x+F(x).$

Renegar’s method was based on a type of interior point methods known as the
_path following method_ , which incrementally minimizing $f_{t}(x)$. Once we
get the minimizer $x_{t}^{*}$ of $f_{t}(x)$ for some large $t$, it is not far
from being optimal, e.g., $c^{\top}x_{t}^{*}\leq\mathrm{OPT}+m/t$.

From Renegar’s framework, the potential function $\Phi_{t}(x):=\|\nabla
f_{t}(x)\|_{H(x)}$ is used to measure how close from $x$ to the minimizer
$x_{t}^{*}$, where $H(x):=\nabla^{2}F(x)$ is the Hessian matrix. $\Phi_{t}(x)$
is always greater or equal than $0$. On the other hand, $\Phi_{t}(x)=0$ means
$x=x_{t}^{*}$. In each iteration, our goal is to increase $t$ to
$t^{\operatorname{new}}$, move $x$ to $x^{\operatorname{new}}$, while keeping
$\Phi_{t^{\operatorname{new}}}(x^{\operatorname{new}})$ small:

$\displaystyle\Phi_{t^{\operatorname{new}}}(x^{\operatorname{new}})=\Phi_{t}(x)+\underbrace{\Phi_{t}(x^{\operatorname{new}})-\Phi_{t}(x)}_{x\mathrm{\leavevmode\nobreak\
move}}+\underbrace{\Phi_{t^{\operatorname{new}}}(x^{\operatorname{new}})-\Phi_{t}(x^{\operatorname{new}})}_{t\mathrm{\leavevmode\nobreak\
move}}.$

We are able to prove the following (see Section 7):

  1. 1.

In the $x$ moving part, let $x^{\operatorname{new}}$ be applying one newton
step from $x$ on function $f_{t}(x)$. Implementing Newton step requires an SDD
solver. After that, we can guarantee the standard quadratic convergence rate
of Newton method:

$\displaystyle\Phi_{t}(x^{\operatorname{new}})\leq\Phi_{t}(x)^{2}.$

  2. 2.

In the $t$ moving part, we multiplicative increase $t$ to
$t^{\operatorname{new}}$. We can guarantee that

$\displaystyle\Phi_{t^{\operatorname{new}}}(x^{\operatorname{new}})\leq\Phi_{t}(x^{\operatorname{new}})+|t^{\operatorname{new}}/t-1|\cdot\sqrt{m}.$

Combining the above two parts, we can set
$t^{\operatorname{new}}=(1+O(1/\sqrt{m}))\cdot t$ so that an
$\widetilde{O}(\sqrt{m})$-iteration IPM follows.

###  2.2 Our implementations in the streaming model

##### Implement IPM in the streaming model.

All the existing IPMs are working in classical model where there is no space
limitation. One of our contribution is to implement IPM in streaming model.

We show that the iterating procedure in IPM can be implemented in the
streaming model and each iteration can be done in $\widetilde{O}(1)$ passes
(see Section 8). Recall that each iteration is performing a Newton step
$(\nabla^{2}F(x))^{-1}\nabla f_{t}(x)\in\mathbb{R}^{n}.$ Note that
$\nabla^{2}F(x)\in\mathbb{R}^{n\times n}$ is an SDD matrix in our setting, and
$\nabla^{2}F(x)$, $\nabla f_{t}(x)$ both rely on the slack variable
$s(x):=Ax-b\in\mathbb{R}^{m}$, but we cannot even store length-$m$ vector $b$
in the space. We address this problem by

  * •

implementing an SDD solver in the streaming model;

  * •

implementing the algorithm in the isolation lemma in the streaming model.

##### Implement the Laplacian solver in the streaming model.

For simplicity, we outline a streaming Laplacian solver here, which gives an
SDD solver by standard reductions (see Section 9). Given a graph $G$, let
$L_{G}$ be the Laplacian matrix of $G$. Let $A$ be the signed edge-vertex
incident matrix, and let $a_{1}^{\top},\cdots,a_{m}^{\top}$ be rows of
$A\in\mathbb{R}^{m\times n}$. Let $s\in\mathbb{R}^{m}$ be slack variable.

Note that $L_{G}$ has $\operatorname{nnz}(L_{G})=O(m)$ so it cannot be stored
in memory, but it can be written as 888We want to remark that another
equivalent definition is $L_{G}=D_{G}-A_{G}$, where
$D_{G}\in\mathbb{R}^{n\times n}$ is a diagonal matrix ($(D_{G})_{i,i}$ is the
weighted degree of vertex $i$) and $A_{G}$ is the adjacency matrix of graph
$G$, where $(A_{G})_{i,j}=(A_{G})_{j,i}$ is the weight of edge $(i,j)$ in the
graph.

$\displaystyle L_{G}=\sum_{i\in[m]}\frac{a_{i}a_{i}^{\top}}{s_{i}^{2}}$

Suppose we have an oracle that outputs $s_{i}$ for any given $i\in[m]$, then
one can do simple things like multiplying $L_{G}$ by a vector $v$, e.g.,
$L_{G}\cdot v$. This can be done since we can read one pass of all edges, and
note that for edge $e_{i}=(u,v)$, we have
$a_{i}=\mathbf{1}_{u}-\mathbf{1}_{v}$ so that
$\frac{a_{i}a_{i}^{\top}}{s_{i}^{2}}\cdot y\in\mathbb{R}^{n}$ can be computed
and stored. In fact, we show that this simple operation suffices: by [KLM+17],
we can compute a sparsifier $H$ of $G$. Using $L_{H}$, the Laplacian of $H$,
as a preconditioner of $L_{G}$, we can do iterative refinement on the
solution:

$\displaystyle x^{(t+1)}:=x^{(t)}+L_{H}^{-1}\cdot(b-L_{G}\cdot x^{(t)}).$

After $\widetilde{O}(1)$ iterations, $x$ will converge to the solution
$x^{*}=L_{G}^{-1}y$ with desired precision.

##### Implement the isolation lemma in the streaming model.

We start with the definition of the isolation lemma (see Section 4 for
details).

######  Definition 2.1 (Isolation lemma).

Given a set system $(S,\mathcal{F})$ where
$\mathcal{F}\subseteq\\{0,1\\}^{S}$. Given weight $w_{i}$ to each element $i$
in $S$, the weight of a set $F$ in $\mathcal{F}$ is defined as $\sum_{i\in
F}w_{i}$. The goal of the isolation lemma is to design a scheme that can
assign weight oblivious to $\mathcal{F}$, such that there is a unique set in
$\mathcal{F}$ that has the minimum (maximum) weight under this assignment.

The isolation lemma always involves randomness since its output is oblivious
to $\mathcal{F}$. The isolation lemma says that if we randomly choose this
weight, then with a good probability the uniqueness is ensured. However, this
does not apply to the streaming setting since this weight vector is of length
$m$, which cannot be stored in memory. The ideal is to _consistently_ output
the same weight vector in different passes.

[CRS95] reduces the number of random bits from $|S|$ to
$\log(|\mathcal{F}|)\leq|S|$. In the case of maximum weight bipartite
matching, we have $|\mathcal{F}|\leq(n+1)^{n}$ since each of the $n$ vertices
can choose none or one of the vertices in the matching. Our streaming
implementation (See Algorithm 2) is to store these
$\log((n+1)^{n})=\widetilde{O}(n)$ random bits, so that we can consistently
output the same weight vector in different passes.

###  2.3 Turn dual to primal

Turning optimal dual solution to optimal primal solution for general LP
requires at least solving a linear system, which takes $O(n^{\omega})$ time
and $O(nm)$ space (Lee, Sidford and Wong [LSW15]), which is unknown to be
implemented in the semi-streaming model even for bipartite matching LP.999In
general, the inverse of a sparse matrix can be dense, which means the standard
Gaussian elimination method for linear system solving does not apply in the
semi-streaming model. We bypass this issue by using the complementary
slackness theorem to highlight $n$ tight dual constraints and therefore
sparsify the original graph from $m$ edges to $n$ edges without losing the
optimal matching. The isolation lemma is the key to implement this idea. To
see why, we give a simple example here showing that the number of highlighted
tight constraints can be as large as $m$ without using the isolation lemma.

Figure 1: The red circle is a minimum vertex cover, which is an optimal dual
solution. The blue edge is a maximum matching, which is an optimal primal
solution. In both examples, the primal and dual satisfy complementary
slackness Eq. (1).

Suppose the graph has a (maximum weight) perfect matching (see Figure 1 for
example). Then the following trivial solution is optimal to the dual LP:
choosing all vertices in $V_{L}$. Let us show what happens when applying the
complementary slackness theorem. The complementary slackness theorem says that
when $y$ is a feasible primal solution and $x$ is a feasible dual solution,
then $y$ is optimal to the primal and $x$ is optimal to the dual if and only
if

$\langle y,Ax-\mathbf{1}_{m}\rangle=0\mathrm{\leavevmode\nobreak\
\leavevmode\nobreak\ and\leavevmode\nobreak\ \leavevmode\nobreak\ }\langle
x,\mathbf{1}_{n}-A^{\top}y\rangle=0.$ (1)

From the above case, we have $Ax-\mathbf{1}_{m}=0$, so the first equality
$\langle y,Ax-\mathbf{1}_{m}\rangle=0$ puts no constraint on $y$. Therefore,
any solution $y\geq\mathbf{0}_{m}$ to the linear system
$a_{i}^{\top}y_{i}=1,\leavevmode\nobreak\ \forall i\in V_{L}$ is an optimal
solution, where $a_{i}$ is the $i$-th column of $A$. Note that this linear
system has $m$ variables and $|V_{L}|$ equations, which is still hard to find
a solution in $\widetilde{O}(n)$ space.9

Roughly speaking, the complementary slackness theorem does not help in this
case since choosing all vertices in $V_{L}$ is always a feasible solution, and
is actually an optimal solution when a perfect matching exists.

Now consider perturbing the primal objective function by some vector
$b\in\mathbb{R}^{m}$ such that the optimal solution to the following primal LP
is unique:

$\displaystyle\textbf{Primal}\leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \max_{y\in\mathbb{R}^{m}}$ $\displaystyle
b^{\top}y,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \text{s.t.
}A^{\top}y\leq\mathbf{1}_{n}\text{ and }y\geq{\bf 0}_{m}.$

Suppose we find the optimal solution $x$ in the dual LP and we want to recover
the optimal solution $y$ in the primal LP. Again by plugging in the
complementary slackness theorem, we get at most $n$ equations from the second
part $\langle x,\mathbf{1}_{n}-A^{\top}y\rangle=0$. Since the optimal $y$ is
unique and $y$ has dimension $m$, the first part $\langle
y,Ax-\mathbf{1}_{m}\rangle$ must contribute at least $m-n$ equations. Note
that these equations have the form

$y_{i}=0,\leavevmode\nobreak\ \forall i\in[m]\text{\leavevmode\nobreak\
s.t.\leavevmode\nobreak\ }(Ax)_{i}-1>0.$

This means that the corresponding edges are _unnecessary_ in order to get one
maximum matching. As a result, we can reduce the number of edges from $m$ to
$n$, then compute a maximum matching in $\widetilde{O}(n)$ space without
reading the stream.

###  2.4 Discussions

We note that going from $\sqrt{m}$ passes to $\sqrt{n}$ passes requires
solving fundamental problems in optimization, (e.g., implementing the Lewis
weights sampling in space-limited models). Given that we already knew that
reachability [LJS19] and single source shortest path [FKM+09, CFCHT20] can be
both solved in $n^{1/2+o(1)}$ passes in semi-streaming, we believe that it is
an interesting open problem to improve the number of passes from $\sqrt{m}$ to
$\sqrt{n}$.

###### Contents

  1. 1 Introduction
    1. 1.1 Our contribution
    2. 1.2 Related work
    3. 1.3 Previous techniques
      1. 1.3.1 Approximation algorithms
      2. 1.3.2 From approximate to exact maximum matching
      3. 1.3.3 A folklore algorithm with $O(n\log n)$ passes
      4. 1.3.4 Other algorithms with more space
  2. 2 Our techniques
    1. 2.1 IPM analysis
    2. 2.2 Our implementations in the streaming model
    3. 2.3 Turn dual to primal
    4. 2.4 Discussions
  3. 3 Notations
  4. 4 Isolation lemma in the streaming model
    1. 4.1 Isolation lemma
    2. 4.2 Algorithms
    3. 4.3 Proof of uniqueness: step 1
    4. 4.4 Proof of uniqueness: step 2
    5. 4.5 Proof of uniqueness: step 3
    6. 4.6 Proof of uniqueness: step 4
  5. 5 Preliminary for IPM
    1. 5.1 Definitions
    2. 5.2 Approximation tools: self-concordant function, near optimal solution and perturbation
  6. 6 Algorithm
  7. 7 Error analysis of IPM
    1. 7.1 Assumptions on parameters
    2. 7.2 Bounding potential function $\Phi$
    3. 7.3 Bounding the movement of $t$
    4. 7.4 Upper bounding the potential function
    5. 7.5 Move both: final
    6. 7.6 Move both: part 1
    7. 7.7 Move both: part 2
    8. 7.8 Move both: part 3
  8. 8 Pass complexity of IPM
  9. 9 SDD solver in the streaming model
    1. 9.1 SDD and Laplacian systems
    2. 9.2 The preconditioner
    3. 9.3 An iterative solver
    4. 9.4 An iterative solver: the space
    5. 9.5 An iterative solver: the accuracy
    6. 9.6 Main result
  10. 10 Minimum vertex cover
    1. 10.1 Our algorithm
    2. 10.2 Correctness of Algorithm 5
    3. 10.3 Pass complexity of Algorithm 5
    4. 10.4 Building blocks
    5. 10.5 A minimum vertex cover solver
  11. 11 Combine
    1. 11.1 Algorithms
    2. 11.2 Pass complexity
    3. 11.3 Correctness
    4. 11.4 Primal to dual
    5. 11.5 Properties of primal and dual LP solutions
  12. A Solver reductions
    1. A.1 From SDDM0 solver to SDD0 solver
    2. A.2 From SDDM solver to SDDM0 solver

##### Roadmap.

The rest of the paper is organized as follows. In Section 3, we define the
basic notations in this paper. In Section 4, we present the generalized
isolation lemma in the semi-streaming model, which will be used to recover the
maximum matching from a minimum vertex cover.

In Section 5, we give some preliminaries for interior point method. In Section
6, we present our streaming algorithm of interior point method. In Section 7,
we present the error analysis of interior point method. In section 8, we show
the pass complexity of our interior point method.

In Section 9, we give an SDD solver in the streaming model, which is a
necessary component of our interior point method. In Section 10, we discuss
the streaming algorithm for minimum vertex cover. In Section 11, we combine
all the pieces together to get the final algorithm for maximum weight
bipartite matching. In Section A, we complement Section 9 by providing two
reductions from weaker solvers to our final SDD0 solver.

##  3 Notations

##### Standard notations.

For a positive integer $n$, we denote $[n]=\\{1,2,\cdots,n\\}$.

We use $\operatorname*{{\mathbb{E}}}[]$ for expectation and $\Pr[]$ for
probability.

For a positive integer $n$, we use $I_{n}$ to denote the identity matrix of
size $n\times n$.

For a vector $v\in\mathbb{R}^{n}$, we use the standard definition of
$\ell_{p}$ norms: $\forall p\geq 1$,
$\|v\|_{p}=(\sum_{i=1}^{n}|v_{i}|^{p})^{1/p}$. Specially,
$\|v\|_{\infty}=\max_{i\in[n]}|v_{i}|$. We use $\|v\|_{0}$ to denote the
number of nonzero entries in vector $v$. We use $\operatorname{supp}(v)$ to
denote the support of vector $v$.

We use ${\bf 0}_{n}$ to denote a length-$n$ vector where every entry is $0$.
We use ${\bf 0}_{m\times n}$ to denote a $m\times n$ matrix where each entry
is $0$. Similarly, we use the notation ${\bf 1}_{n}$ and ${\bf 1}_{m\times
n}$.

##### Matrix operators.

For a square matrix $A$, we use $\operatorname{tr}[A]$ to denote its trace.
For a square and full rank matrix $A$, we use $A^{-1}$ denote the true inverse
of $A$. For a matrix $A$, we use $A^{\dagger}$ to denote its pseudo inverse.
We say a square matrix $A$ is positive definite, if for all $x$,
$x^{\top}Ax>0$. We say a square matrix $A$ is positive semi-definite, if for
all $x$, $x^{\top}Ax\geq 0$. We use $\succeq$ and $\preceq$ to denote the
p.s.d. ordering. For example, we say $A\succeq B$, if $x^{\top}Ax\geq
x^{\top}Bx,\forall x$.

##### Matrix norms.

For a matrix $A$, we use $\|A\|_{1}$ to denote its entry-wise $\ell_{1}$ norm,
i.e., $\|A\|_{1}=\sum_{i,j}|A_{i,j}|$. We use $\|A\|_{F}$ to denote its
Frobenius norm $\|A\|_{F}=(\sum_{i,j}A_{i,j}^{2})^{1/2}$. We use $\|A\|$ to
denote its spectral/operator norm.

##### Matrix approximation.

Let $A,B\in\mathbb{R}^{n\times n}$ be positive semi-definite matrix. Let
$\epsilon\in(0,1)$. We say $A\approx_{\epsilon}B$ if

$(1-\epsilon)x^{\top}Ax\leq
x^{\top}Bx\leq(1+\epsilon)x^{\top}Ax,\leavevmode\nobreak\ \leavevmode\nobreak\
\forall x\in\mathbb{R}^{n}.$

Note that if we have $A\approx_{\epsilon}B$, then
$(1-\epsilon)\|x\|_{A}\leq\|x\|_{B}\leq(1+\epsilon)\|x\|_{A}$ for all
$x\in\mathbb{R}^{n}$.

##### The SDD matrix-related definitions.

######  Definition 3.1 (Edge-vertex incident matrix).

Let $G=(V_{L},V_{R},E)$ be a connected undirected bipartite graph. The
_(unsigned) edge-vertex incidence matrix_ is denoted as follows

$\displaystyle B(e,v)=\begin{cases}1,&\mathrm{\leavevmode\nobreak\
if\leavevmode\nobreak\ }v\mathrm{\leavevmode\nobreak\
incident\leavevmode\nobreak\ to\leavevmode\nobreak\ }e;\\\
0,&\mathrm{\leavevmode\nobreak\ otherwise}.\end{cases}$

######  Definition 3.2 (signed edge-vertex incident matrix).

Let $G=(V_{L},V_{R},E)$ be a connected directed bipartite graph where all
edges orientate from $V_{L}$ to $V_{R}$. The _signed edge-vertex incidence
matrix_ is denoted as follows

$\displaystyle B(e,v)=\begin{cases}+1,&\mathrm{\leavevmode\nobreak\
if\leavevmode\nobreak\ }v\mathrm{\leavevmode\nobreak\
incident\leavevmode\nobreak\ to\leavevmode\nobreak\
}e\mathrm{\leavevmode\nobreak\ and\leavevmode\nobreak\ }v\in V_{L};\\\
-1,&\mathrm{\leavevmode\nobreak\ if\leavevmode\nobreak\
}v\mathrm{\leavevmode\nobreak\ incident\leavevmode\nobreak\
to\leavevmode\nobreak\ }e\mathrm{\leavevmode\nobreak\ and\leavevmode\nobreak\
}v\in V_{R};\\\ 0,&\mathrm{\leavevmode\nobreak\ otherwise}.\end{cases}$

######  Definition 3.3 (SDDM, SDD matrix).

A square matrix $A$ is _weakly diagonally dominant_ if $A_{i,i}\geq\sum_{j\neq
i}\left|A_{i,j}\right|$ for $i$, and is _strictly diagonally dominant_ if
$A_{i,i}>\sum_{j\neq i}\left|A_{i,j}\right|$ for $i$. A matrix $A$ is SDD0 if
it is symmetric and weakly diagonally dominant, and is SDD if it is symmetric
and strictly diagonally dominant. A matrix $A$ is SDDM0 if it is SDD and
$A_{i,j}\leq 0$ for all $i\neq j$. An SDDM0 matrix is SDDM if it is strictly
diagonally dominant.

######  Fact 3.4.

An SDDM0 matrix must be positive semi-definite. If an SDDM0 matrix has zero
row-sums, then it is a Laplacian matrix. If an SDDM0 matrix has at least one
positive row-sum, then it is an SDDM matrix and positive definite.

##### Bit complexity.

Given a linear programming

$\displaystyle\leavevmode\nobreak\ \min_{x\in\mathbb{R}^{n}}c^{\top}x$ (2)
s.t. $\displaystyle\leavevmode\nobreak\ Ax\geq b,$

where $A\in\operatorname*{\mathbb{Z}}^{m\times n}$,
$b\in\operatorname*{\mathbb{Z}}^{m}$, $c\in\operatorname*{\mathbb{Z}}^{n}$ all
having integer coefficient.

The bit complexity $L$ is defined as

$\displaystyle
L:=\log(m)+\log(1+\mathrm{d}_{\mathrm{max}}(A))+\log(1+\max\\{\|c\|_{\infty},\|b\|_{\infty}\\}),$

where $\mathrm{d}_{\mathrm{max}}(A)$ denotes the largest absolute value of the
determinant of a square sub-matrix of $A$.

It is well known that $\operatorname{poly}(L)$-bit precision is sufficient to
implement an IPM (e.g., see [DS08] and the references therein). This is
because the absolute values of all intermediate arithmetic results are within
$[2^{-\operatorname{poly}(L)},2^{\operatorname{poly}(L)}]$, and the errors in
all the approximations are at least $1/\operatorname{poly}(n)$. Therefore,
truncating all the arithmetic results to $\operatorname{poly}(L)$ bits for
some sufficiently large polynomial preserves all the error parameters and thus
the same analysis holds.

We will need the following tools in our later analysis.

######  Lemma 3.5 ([HT56]).

Let $A\in\mathbb{R}^{m\times n}$ be the unsigned edge-vertex incident matrix
of a bipartite graph $G=(V_{L},V_{R},E)$. Let $S:=\\{Ax\leq b\mid
x\geq\mathbf{0}_{m}\\}$, $S^{\prime}:=\\{A^{\top}y\leq b^{\prime}\mid
y\geq\mathbf{0}_{n}\\}$, where $b\in\operatorname*{\mathbb{Z}}^{m}$,
$b^{\prime}\in\operatorname*{\mathbb{Z}}^{n}$. Then both $A$ and $A^{\top}$
are totally unimodular, i.e., all square submatrices of them have determinants
of $0,1,-1$, and all extreme points of $S$ and $S^{\prime}$ are integral.

##  4 Isolation lemma in the streaming model

[CRS95] already shows how to implement isolation lemma using a small amount of
randomness, and in particular, in our application, the amount of randomness is
$O(n)$ and therefore just fits into our memory. However, since their focus is
on the number of randomness, they still use extra space that we cannot
tolerate. In this section, we make their algorithm into an oracle. This oracle
only stores the random seed and performs exactly the same as the original
algorithm, so that we can use it in the streaming model. Formally, our result
is stated in Lemma 4.1.

Before we prove, we give a simple notation here.

For a vector $w\in\mathbb{Z}^{n}$ and a set $S\subseteq[n]$, we denote
$w_{S}:=\sum_{i\in S}w_{i}$. We can define $w_{S}:=\sum_{i\in S}w_{i}$
similarly when $w_{i}:\mathbb{Z}^{t}\rightarrow\mathbb{Z}$ is a function.

######  Lemma 4.1 (Streaming implementation of isolation lemma).

Let $n,\mathcal{F},Z,w$ be define as in Lemma 4.2. The Algorithm 1 in Lemma
4.2 can be implemented into such an oracle $\mathcal{I}$: $\mathcal{I}$ can
output $w_{i}$ given any $i\in[n]$. Furthermore, $\mathcal{I}$ uses
$O(\log(Z)+\log(n))$ space.

###### Proof.

Our $\mathcal{I}$ is the streaming implementation of Algorithm 2. It is easy
to see that after running the procedure Initialize, the procedure
$\textsc{Query}(i)$ will output $w_{i}$ given any $i\in[n]$.

This oracle stores $m,Z\in\mathbb{N}$ and $r\in\mathbb{N}^{t}$ in memory. Note
that $m,Z$ can be stored in $O(\log(Z)+\log(n))$ bits, and $r$ can be stored
in

$O(t\cdot\log
n)=O(\lceil\log(m)/\log(n)\rceil\cdot\log(n))=O(\log(Z)+\log(n))$

bits. And it is easy to see in Query, all computation can be done within
$O(\log(Z)+\log(n))$ space. ∎

The rest of this section is organized as follows: In Section 4.1, we the
generalized isolation lemma. In Section 4.2, we provide our streaming
algorithm (in fact a data-structure). In Section 4.3, Section 4.4, Section 4.5
and Section 4.6, we provide the proof details of the uniqueness.

###  4.1 Isolation lemma

We state the generalized isolation lemma from previous work [CRS95].

######  Lemma 4.2 (Generalized isolation lemma [CRS95]).

Fix $n\in\mathbb{N}$. Fix an unknown family $\mathcal{F}\subseteq[n]$. Let $Z$
denote a positive integer such that $Z\geq|\mathcal{F}|$, there exists an
algorithm (Algorithm 1) that uses $O(\log(Z)+\log(n))$ random bits to output a
vector $w\in[0,n^{7}]^{n}$, such that with probability at least $1/4$, there
is a unique set $S^{*}\in\mathcal{F}$ that has minimum weight $w_{S^{*}}$.

###### Proof.

The proof is already done in [CRS95]. For the completeness, we rewrite their
proof here. By Lemma 4.4, with probability at least $1/2$, all sets
$S\in\mathcal{F}$ has distinct $w^{(2)}_{S}$. Under the event of we success,
by Lemma 4.5, we get that all sets $S\in\mathcal{F}$ has distinct
$w^{(3)}_{S}$. Then by Lemma 4.6, we get that with probability at least $1/2$,
the output $w$ will give a unique minimum set in $\mathcal{F}$. Since the two
events of success are independent, the final success probability is at least
$1/4$. ∎

###  4.2 Algorithms

We present an algorithm that implements the isolation lemma in the streaming
model. We give the original implementation appeared in [CRS95] in Algorithm 1.
Then in Algorithm 2, we show how to implement the algorithm in the streaming
model.

We explain the Algorithm 2. Algorithm 2 is in fact a data-structure and has
three parts. The first part is all the members. The second part is a function
for initialization (see Initialize). The third part is a function for query
(see Query). The function Initialize is initializing variables $m,t$ and
vector $r$. The Query function takes $i\in[n]$ as input and output an integer
$w^{(4)}$. We want to remark that this data-structure is a static data-
structure. Therefore, it doesn’t need an Update function.

Algorithm 1 Conceptual implementation of the isolation lemma. Algorithm 2 is
the implementation of this algorithm in streaming model.

1:procedure Isolation($n,Z\in\mathbb{N}$) $\triangleright$
$Z\geq|\mathcal{F}|$, Lemma 4.2

2: $w_{i}^{(1)}\leftarrow 2^{i}$, $\forall i\in[n]$ $\triangleright$
$w^{(1)}\in\mathbb{N}^{n}$

3: Choose $m$ uniformly at random from $\\{1,2,\cdots,(2nZ^{2})^{2}\\}$.
$\triangleright$ $m\leq 4n^{2}Z^{4}$

4: For each $i$, define $w^{(2)}_{i}\leftarrow w^{(1)}_{i}\mod m$.
$\triangleright$ $w^{(2)}\in[m]^{n}$

5: $t\leftarrow\lceil\log(m)/\log(n)\rceil$

6: for $i=1\to n$ do

7: $\overline{b_{i,t-1},\cdots,b_{i,1},b_{i,0}}\leftarrow w_{i}^{(2)}$
$\triangleright$ Write $w_{i}^{(2)}$ in base $n$. $b_{i,j}\in[n]$ are digits.

8: $\triangleright$ Note that $t$ is an upper bound on the length

9: $w^{(3)}_{i}(y_{0},\cdots,y_{t-1})\leftarrow\sum_{j=0}^{t-1}b_{i,j}\cdot
y_{j}$ $\triangleright$ $w^{(3)}_{i}:\mathbb{Z}^{t}\rightarrow\mathbb{Z}$ is a
linear form.

10: end for

11: Choose $r_{0},\cdots,r_{t-1}$ uniformly and independently at random from
$\\{1,2,\cdots,n^{5}\\}$.

12: $w^{(4)}_{i}\leftarrow w^{(3)}_{i}(r_{0},\cdots,r_{t-1})$, $\forall
i\in[n]$

13: return $w^{(4)}\in\mathbb{N}^{n}$.

14:end procedure

Algorithm 2 Data structure: Isolation(). This algorithm is the streaming
implementation of Algorithm 1.

1:data structure $\triangleright$ Lemma 4.1

2:members

3: $n,t\in\mathbb{N}$ $\triangleright$ $t=O(\log(Z)/\log(n))$

4: $m,Z\in\mathbb{N}$ $\triangleright$ $m=O(n^{2}Z^{4})$, $Z\geq|\mathcal{F}|$

5: $r_{0},\cdots,r_{t-1}\in\mathbb{N}$ $\triangleright$ $r_{i}\leq n^{5}$

6:end members

7:procedure Initialize($n,Z\in\mathbb{N}$) $\triangleright$ Initialization

8: Choose $m\in\mathbb{N}$ uniformly at random from
$\\{1,2,\cdots,(2nZ^{2})^{2}\\}$

9: $t\leftarrow\lceil\log m/\log n\rceil$ $\triangleright$ $t\in\mathbb{N}$

10: Choose $r_{0},\cdots,r_{t-1}$ uniformly and independently at random from
$\\{1,2,\cdots,n^{5}\\}$ $\triangleright$ $r\in\mathbb{N}^{t}$

11:end procedure

12:procedure Query($i\in[n]$)

13: $w^{(1)}\leftarrow 2^{i}$ $\triangleright$ $w^{(1)}\in\mathbb{N}$,
$w^{(1)}\leq 2^{n}$

14: $w^{(2)}\leftarrow w^{(1)}_{i}\mod m$ $\triangleright$
$w^{(2)}\in\mathbb{N}$, $w^{(2)}\leq 2^{n}$

15: $\overline{b_{t-1},\cdots,b_{1},b_{0}}\leftarrow w^{(2)}$ $\triangleright$
$b_{j}\in[n]$, $\forall j\in[t]$

16: $\triangleright$ Write $w^{(2)}$ in base $n$. Note that $t$ is an upper
bound on the length

17: $w^{(4)}\leftarrow\sum_{j=0}^{t-1}b_{j}\cdot r_{j}$

18: return $w^{(4)}$ $\triangleright$ $w^{(4)}\leq n^{7}$

19:end procedure

20:end data structure

###  4.3 Proof of uniqueness: step 1

The goal of this section is to prove Lemma 4.3. The following lemma is from
[Thr93].

######  Lemma 4.3 (Step 1, [Thr93]).

Let $L\geq 100$ and let $S$ be any subset of $\\{1,\cdots,L^{2}\\}$ such that
$|S|\geq\frac{1}{2}L^{2}$. Then, the least common multiple of the elements in
$S$ exceeds $2^{L}$.

###  4.4 Proof of uniqueness: step 2

The goal of this section is to prove Lemma 4.4.

######  Lemma 4.4 (Step 2).

With probability at least $1/2$, all sets $S\in\mathcal{F}$ have distinct
weights $w^{(2)}_{S}$.

###### Proof.

We write $\mathcal{F}=\\{S_{1},\cdots,S_{k}\\}$ where $k\leq Z$. Define

$I:=\prod_{1\leq i<j\leq k}\big{|}w^{(1)}_{S_{i}}-w^{(1)}_{S_{j}}\big{|}.$

First note that $I\neq 0$ since every set $S\in\mathcal{F}$ have distinct
weights $w^{(1)}_{S}$ by definition of $w^{(1)}_{i}=2^{i}$ (Line 2, Algorithm
1).

Next, we give an upper bound on $I$. For each pair of $1\leq i<j\leq k$,
$|w^{(1)}_{S_{i}}-w^{(1)}_{S_{j}}|\leq 2^{n+1}$. There are totally $Z^{2}$
pairs, so $I<2^{2nZ^{2}}$.

Let $L=2nZ^{2}$. Let $S=\\{1,\cdots,L^{2}\\}$ be all the possible choices of
$m$. We have that at least half choices of $m$ satisfies $I\mod m\neq 0$,
since otherwise by applying Lemma 4.3, $I$ is at least the least common
multiplier of half numbers of $S$, i.e., $I\geq 2^{L}$, contradicting with the
upper bound on $I<2^{2nZ^{2}}$.

Therefore, with probability at least $1/2$ (the randomness is over the choice
of $m$), $I\mod m\neq 0$, which means that all sets $S\in\mathcal{F}$ have
distinct weight $w^{(2)}_{S}$ by our definition of
$w^{(2)}_{i}=w^{(1)}_{i}\mod m$ (Line 4, Algorithm 1). ∎

###  4.5 Proof of uniqueness: step 3

The goal of this section is to prove Lemma 4.5.

######  Lemma 4.5 (Step 3).

If $w^{(2)}_{S}$ are distinct for all $S\in\mathcal{F}$, then the linear form
$w^{(3)}_{S}$ are all distinct for all $S\in\mathcal{F}$.

###### Proof.

Use proof by contradiction. Suppose there exists two distinct sets
$S_{1},S_{2}\in\mathcal{F}$ that $w^{(3)}_{S_{1}}=w^{(3)}_{S_{2}}$. Let
$y_{i}=n^{i}$, $\forall i\in[t]$. We will have

$w^{(2)}_{S_{1}}=w^{(3)}_{S_{1}}(y_{0},\cdots,y_{t-1})=w^{(3)}_{S_{2}}(y_{0},\cdots,y_{t-1})=w^{(2)}_{S_{2}}$

by definition of $w^{(3)}_{i}$ (Line 9), which is a contradict to our
assumption that all $w^{(2)}_{S}$ are distinct. ∎

###  4.6 Proof of uniqueness: step 4

The goal of this section is to prove Lemma 4.6.

######  Lemma 4.6 (Step 4).

Let $\mathcal{C}$ be any collection of distinct linear forms over at most $t$
variables $y=y_{0},\cdots,y_{t-1}$ with coefficients in
$\\{0,1,\cdots,n-1\\}$. Choose a random $r=r_{0},\cdots,r_{t-1}$ by assigning
each $r_{i}$ uniformly and independently from $[n^{2}\cdot t]$. Then in the
assignment $y=r$ there will be a unique linear form with minimum value, with
probability at least $1/2$.

###### Proof.

We call a variable $y_{i}$ to be singular under an assignment
$r\in\operatorname*{\mathbb{Z}}^{t}$ if there exists two minimum linear forms
in $\mathcal{C}$ under assignment $r$. Then an assignment $r$ gives unique
minimum linear form if and only if no variable $y_{i}$ is singular under this
assignment. We will calculate the probability of $y_{i}$ being singular under
random assignment $r$ and then take union bound over every $y_{i}$.

For each $y_{i}$, fix all
$r_{0},\cdots,r_{i-1},r_{i+1},\cdots,r_{t-1}=a_{0},\cdots,a_{i-1},a_{i+1},\cdots,a_{t-1}$
other than $r_{i}$. Now, every linear form $f$ under this partial assignment
can be written as $a_{f}+b_{f}\cdot y_{i}$ with $b_{f}<n$. We split
$\mathcal{C}$ into $n$ classes $\mathcal{C}_{0},\cdots,\mathcal{C}_{n-1}$
where $\mathcal{C}_{j}$ contains all linear forms with $b_{f}=j$. Let $p_{j}$
be the minimum $a_{f}$ among all linear forms in $\mathcal{C}_{j}$. According
to the definition of singular, $y_{i}$ is singular on assignment $r_{i}$ if
and only if the minimum value in the list

$\\{p_{0},p_{1}+r_{i},p_{2}+2r_{i},\cdots,p_{n-1}+(n-1)\cdot r_{i}\\}$

is not unique, which is upper bounded by the probability that the elements in
the list has a collision. Since every pair of elements in the list can have at
most one choice of $r_{i}$ such that they are equal, we have

$\displaystyle\Pr_{r_{i}}[y_{i}\text{\leavevmode\nobreak\ is singular}\mid
r_{j}=a_{j},\leavevmode\nobreak\ \forall j\neq i]\leq$
$\displaystyle\leavevmode\nobreak\ \Pr_{r_{i}}[\exists
l,m\text{\leavevmode\nobreak\ s.t.\leavevmode\nobreak\ }p_{l}+l\cdot
r_{i}=p_{m}+m\cdot r_{i}]$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \sum_{l\neq m}\Pr_{r_{i}}[p_{l}+l\cdot
r_{i}=p_{m}+m\cdot r_{i}]$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \binom{n}{2}\cdot\frac{1}{n^{2}t}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\ \frac{1}{2t}.$

Finally, by a union bound on the events that every $y_{i}$ is not singular,
the conclusion follows with probability at least $1-t\cdot\frac{1}{2t}=1/2$. ∎

##  5 Preliminary for IPM

Since IPM was proposed by Karmarkar [Kar84] in 1984, it becomes a very popular
method for analyzing the running time of linear programming or linear
programming type problems [Vai87, Ren88, Vai89, NN89, NN92, NN94, Ans00,
Ren01, DS08, Mad13, LS13, LS14, LS15, Mad16, CLS19, LSZ19, BLSS20, Bra20,
JKL+20, LS20a, BLN+20, LS20b, Ye20, HJST21, SY21, JSWZ21, DLY21]. In this
section, we will focus on introducing some mathematical background of IPM for
LP. Let us consider the linear programming

$\displaystyle\min_{x\in\mathbb{R}^{n}}c^{\top}x,$
$\displaystyle\mathrm{\leavevmode\nobreak\ \leavevmode\nobreak\
s.t.}\leavevmode\nobreak\ Ax\geq b$

where $A\in\mathbb{R}^{m\times n}$, $b\in\mathbb{R}^{m}$,
$c\in\mathbb{R}^{n}$. In the rest of this section, all the notations and
discussions will be based on the above LP formulation. For convenient, we
denote $a_{1}^{\top},\cdots,a_{m}^{\top}$ as the row vectors of matrix
$A\in\mathbb{R}^{m\times n}$.

This section has two parts. In Section 5.1, we provide several definitions
related to IPM analysis. In Section 5.2, we state several tools from previous
work.

###  5.1 Definitions

In this section, we give definitions for barrier function, feasible solution,
slack variables, gradient, hessian, perturbed objective function, Netwon step,
and potential functions.

Statement | Comments
---|---
Definition 5.1 | Barrier function
Definition 5.2 | Feasible solution
Definition 5.3 | Slack variables
Definition 5.4 | Gradient and Hessian with respect to barrier function
Definition 5.5 | Perturbed objective function
Definition 5.6 | Newton step
Definition 5.7 |  $\Phi_{t}$, depends on $t$
Definition 5.8 |  $\Psi$, does not depend on $t$
Table 1: IPM related definitions.

We define barrier function:

######  Definition 5.1 (Barrier function).

Define the _logarithmic barrier function_
$F(x):\mathbb{R}^{n}\rightarrow\mathbb{R}$ as follows:

$\displaystyle F(x)\coloneqq-\sum_{i\in[m]}\ln(a_{i}^{\top}x-b_{i}).$ (3)

We define feasible solution as follows:

######  Definition 5.2 (Feasible solution).

For any $x\in\mathbb{R}^{n}$, we say $x$ is feasible if for all $i\in[m]$,
$a_{i}^{\top}x>b_{i}$.

We define slack variables:

######  Definition 5.3 (Slack).

For simplicity, we also define the _slack_ $s_{i}(x)\coloneqq
a_{i}^{\top}x-b_{i}\in\mathbb{R}$ for all $x\in\mathbb{R}^{n}$ and $i\in[m]$.

######  Definition 5.4 (Gradient and Hessian of barrier function).

We define the gradient $g(x):\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}$ and
Hessian $H(x):\mathbb{R}^{n}\rightarrow\mathbb{R}^{n\times n}$ of the barrier
function $F(x)$ as follows:

$\displaystyle g(x):=$ $\displaystyle\leavevmode\nobreak\ \nabla
F(x)=-\sum_{i\in[m]}\frac{a_{i}}{s_{i}(x)}\in\mathbb{R}^{n}$ (4)
$\displaystyle H(x):=$ $\displaystyle\leavevmode\nobreak\
\nabla^{2}F(x)=\sum_{i\in[m]}\frac{a_{i}a_{i}^{\top}}{s_{i}(x)^{2}}\in\mathbb{R}^{n\times
n}.$ (5)

######  Definition 5.5 (Perturbed objective function).

Define the _perturbed objective function_
$f_{t}:\mathbb{R}^{n}\rightarrow\mathbb{R}$, where $t>0$ is a parameter:

$\displaystyle f_{t}(x):=t\cdot c^{\top}x+F(x).$ (6)

######  Definition 5.6 (Newton step).

Given a scalar $t>0$ and a vector $c\in\mathbb{R}^{n}$, define the _Newton
step_ at a feasible point $x\in\mathbb{R}^{n}$ to be

$\displaystyle\delta_{x}:=-(H(x))^{{\dagger}}\cdot\nabla
f_{t}(x)=-(H(x))^{{\dagger}}\cdot(t\cdot c+\nabla
F(x))=-(H(x))^{{\dagger}}\cdot(t\cdot c+g(x)).$ (7)

We define our potential function $\Phi$.

######  Definition 5.7 (Definition of potential $\Phi$).

Given $t>0$ and feasible $x,y\in\mathbb{R}$, we define function
$\Phi_{t}:\mathbb{R}^{n}\times\mathbb{R}^{n}\rightarrow\mathbb{R}$:

$\displaystyle\Phi_{t}(x,y):=\|\nabla f_{t}(x)\|_{H(y)^{-1}}.$

In addition, we also need to define a potential function $\Psi$.

######  Definition 5.8 (Helper of potential function $\Psi$).

Given feasible $x,y\in\mathbb{R}$, we define function
$\Psi:\mathbb{R}^{n}\times\mathbb{R}^{n}\rightarrow\mathbb{R}$:

$\displaystyle\Psi(x,y):=\|g(x)\|_{H(y)^{-1}}.$

###  5.2 Approximation tools: self-concordant function, near optimal solution
and perturbation

We list three types of tools from previous literature. The first one is tools
for analyzing property of self-concordant function. The second one is about
the property of near optimal solution. We remark that those ideas was firstly
proposed in [Ren88], we cite the more cleaner statement from Regenar’s book
[Ren01]. The last one is about perturbation.

######  Lemma 5.9 (Hessian approximation, [Ren01]).

Let $f$ be a self-concordant function with domain $D_{f}$. Define
$H(x):=\nabla^{2}f(x)$. For any two feasible point $x,z\in D_{f}$, if
$\|x-z\|_{H(x)}\leq 1/4$, then we have

$(1-\|x-z\|_{H(x)})^{2}\cdot H(z)\preceq
H(x)\preceq(1-\|x-z\|_{H(x)})^{-2}\cdot H(z).$

######  Lemma 5.10 ([Ren01]).

Let $f$ be a self-concordant function with domain $D_{f}$. For some $x\in
D_{f}$, define $H(x):=\nabla^{2}f(x)$ and $g(x):=\nabla f(x)$. Define
$n(x):=H(x)^{-1}g(x)$ as the Newton step at $x$. If $\|n(x)\|_{H(x)}\leq 1/4$,
then $f$ has a minimizer $z$ that

$\|z-x\|_{H(x)}\leq\|n(x)\|_{H(x)}+\frac{3\|n(x)\|_{H(x)}^{2}}{(1-\|n(x)\|_{H(x)})^{3}}$

######  Lemma 5.11 ([Ren01]).

Let $f$ be a self-concordant function with domain $D_{f}$. Let
$c\in\mathbb{R}^{n}$ be the objective vector. Define $\mathrm{OPT}:=\min_{x\in
D_{f}}\langle c,x\rangle$. Denote $H(x):=\nabla^{2}f(x)$. For any $x,y\in
D_{f}$, we have

$\frac{\langle c,y\rangle-\mathrm{OPT}}{\langle c,x\rangle-\mathrm{OPT}}\leq
1+\|y-x\|_{H(x)}.$

###### Proof.

Define $c_{x}:=H(x)^{-1}c$. Define $B(x)\subseteq\mathbb{R}^{n}$ be the
ellipsoid $B(x):=\\{y\in\mathbb{R}^{n}\mid\|y-x\|_{H(x)}<1\\}$. Since $f$ is
self-concordant, we know $B(x)\subseteq D_{f}$. Therefore for all $0\leq
t<1/\|c_{x}\|_{H(x)}$, $x-t\cdot c_{x}\in D_{f}$. This means

$\langle c,x\rangle-t\cdot\|c_{x}\|_{H(x)}^{2}\geq\mathrm{OPT}.$

By letting $t\rightarrow 1/\|c_{x}\|_{H(x)}$, we get

$\|c_{x}\|_{H(x)}\leq\langle c,x\rangle-\mathrm{OPT}.$

Hence

$\displaystyle\frac{\langle c,y\rangle-\mathrm{OPT}}{\langle
c,x\rangle-\mathrm{OPT}}=$ $\displaystyle\leavevmode\nobreak\ 1+\frac{\langle
c_{x},y-x\rangle_{H(x)}}{\langle c,x\rangle-\mathrm{OPT}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
1+\frac{\|c_{x}\|_{H(x)}\|y-x\|_{H(x)}}{\langle c,x\rangle-\mathrm{OPT}}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\ 1+\|y-x\|_{H(x)},$

where the first step is by the definition of inner product that for some
positive semi-definite matrix $S\in\mathbb{R}^{n\times n}$, $\langle
x,y\rangle_{S}:=\langle Sx,y\rangle$, the second step is by Cauchy-Schwarz
inequality.

∎

######  Lemma 5.12 ($x_{t}^{*}$ is nearly-optimal to $x^{*}$ [Ren01]).

Given linear programming

$\min_{x\in\mathbb{R}^{n},Ax\geq b}c^{\top}x,$

where $A\in\mathbb{R}^{m\times n}$, $b\in\mathbb{R}^{m}$,
$c\in\mathbb{R}^{n}$. Let $x^{*}\in\mathbb{R}^{n}$ be the optimal solution of
the above LP. For every $t>0$, let $x^{*}_{t}\in\mathbb{R}^{n}$ be the $x$
that minimizes $f_{t}(x)$ (Definition 5.5), it must be that

$\displaystyle c^{\top}x^{*}_{t}-c^{\top}x^{*}<m/t.$

###### Proof.

Since $x^{*}_{t}\in\mathbb{R}^{n}$ is the minimizer of $f_{t}(x)$, we have
$\nabla f_{t}(x^{*}_{t})=t\cdot c+g(x_{t}^{*})=\mathbf{0}_{n}$. By
$g(x_{t}^{*})=-t\cdot c\in\mathbb{R}^{n}$, we get

$\displaystyle
c^{\top}x^{*}_{t}-c^{\top}x^{*}=(g(x^{*}_{t}))^{\top}\cdot(x^{*}-x^{*}_{t})/t.$
(8)

Note that $x^{*}$, $x^{*}_{t}\in\mathbb{R}^{n}$ are both feasible, so we have

$\displaystyle(g(x^{*}_{t}))^{\top}\cdot(x^{*}-x^{*}_{t})=$
$\displaystyle\leavevmode\nobreak\
\sum_{i\in[m]}\frac{a_{i}^{\top}(x^{*}-x^{*}_{t})}{s_{i}(x^{*}_{t})}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\sum_{i\in[m]}\frac{s_{i}(x^{*}_{t})-s_{i}(x^{*})}{s_{i}(x^{*}_{t})}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
m-\sum_{i\in[m]}\frac{s_{i}(x^{*})}{s_{i}(x^{*}_{t})}$ $\displaystyle<$
$\displaystyle\leavevmode\nobreak\ m.$

where the first step is by Definition of $g(x)$ (Definition 5.4), the second
step is by Definition of $s_{i}(x)$ (Definition 5.3) that
$s_{i}(x)=a_{i}^{\top}x-b$.

Combining with (Eq. (8)), the lemma follows. ∎

######  Lemma 5.13 (Nearly-optimal output: value version [Ren01]).

Given linear programming

$\min_{x\in\mathbb{R}^{n},Ax\geq b}c^{\top}x,$

where $A\in\mathbb{R}^{m\times n}$, $b\in\mathbb{R}^{m}$,
$c\in\mathbb{R}^{n}$. Let $x^{*}\in\mathbb{R}^{n}$ be the optimal solution of
the above LP. Let $\epsilon_{N}\in(0,1/10)$. If for some $t>0$ and
$x\in\mathbb{R}^{n}$ we have $\Phi_{t}(x,x)\leq\epsilon_{N}$ (Definition 5.7),
then we have

$\displaystyle c^{\top}x-c^{\top}x^{*}\leq\frac{m}{t}\cdot(1+2\epsilon_{N}).$

###### Proof.

Denote $\mathrm{OPT}:=c^{\top}x^{*}$ as the best value we can get in LP.

Let $z$ be the minimizer of $f_{t}$. From $\Phi_{t}(x,x)\leq\epsilon_{N}$, we
know $\|n(x)\|_{H(x)}\leq\epsilon_{N}$, where $n(x):=H(x)^{-1}\nabla
f_{t}(x)$.

By Lemma 5.10, we have

$\displaystyle\|z-x\|_{H(x)}\leq 2\epsilon_{N}.$ (9)

Plug Eq. (9) into Lemma 5.11, we get

$\displaystyle\frac{c^{\top}z-\mathrm{OPT}}{c^{\top}x-\mathrm{OPT}}\leq
1+2\epsilon_{N}.$

By Lemma 5.12, we have

$\displaystyle c^{\top}z\leq\mathrm{OPT}+m/t.$

By combining the above three inequalities, we get

$\displaystyle c^{\top}x-\mathrm{OPT}\leq\frac{m}{t}\cdot(1+2\epsilon_{N}).$

∎

######  Lemma 5.14 (Perturbed linear system [BLSS20]).

Let $H\in\mathbb{R}^{d\times d}$ be a symmetric positive definite matrix. Let
$\epsilon\in(0,1/80)$. Let $b,v\in\mathbb{R}^{d}$ such that
$\|v\|_{H^{-1}}\leq\epsilon\|b\|_{H^{-1}}$. Let
$y=H^{-1}(b+v)\in\mathbb{R}^{d}$. Then there is a symmetric matrix
$\Delta\in\mathbb{R}^{d\times d}$, such that
$y=(H+\Delta)^{-1}b\in\mathbb{R}^{d}$ and

$(1+20\epsilon)^{-1}H\preceq H+\Delta\preceq(1+20\epsilon)H$

##  6 Algorithm

In previous sections, we already define the basic notations (Section 3). In
Section 4, we present the generalized isolation lemma in the semi-streaming
model, which will be used to recover the maximum matching from a minimum
vertex cover. In Section 5, we give some preliminaries for interior point
method.

In this section, we present our IPM algorithm.

Later in Section 7, we will give the error analysis of the output of this
algorithm. In Section 8, we bound the pass number while executing this
algorithm in the streaming model. Specifically, we implement and analyze
streaming SDD solver in Section 9 which is used in Line 16 of this algorithm.

Algorithm 3 Our interior point method

1:procedure
InteriorPointMethod($A,b,c,t_{\mathrm{start}},x_{\mathrm{start}},t_{\mathrm{final}},\epsilon_{\Phi},o$)
$\triangleright$ Lemma 7.1

2: $m,n\leftarrow$ dimensions of $A$

3: $\triangleright$ $A\in\mathbb{R}^{m\times n}$ is the input matrix,
$b\in\mathbb{R}^{m}$, $c\in\mathbb{R}^{n}$

4: $\triangleright$ $t_{\mathrm{start}},x_{\mathrm{start}}$ satisfy
$\Phi_{t_{\mathrm{start}}}(x_{\mathrm{start}},x_{\mathrm{start}})\leq\epsilon_{\Phi}$

5: $\triangleright$ $t_{\mathrm{final}}\in\mathbb{R}$ is the final goal of $t$

6: $\triangleright$ $o\in\\{0,1\\}$. If $o$ is $0$, the algorithm decreases
$t$, otherwise the algorithm increases $t$

7: $\epsilon_{x}\leftarrow 1/100$

8: $\epsilon_{t}\leftarrow 1/(100\sqrt{m})$

9: $\epsilon_{\Phi}\leftarrow\epsilon_{\Phi}$$\triangleright$
$\epsilon_{\Phi}\leq 1/100$ is given from outside

10: $t\leftarrow t_{\mathrm{start}}$

11: $x\leftarrow x_{\mathrm{start}}$

12: $T\leftarrow
O(\sqrt{m})\cdot|\log(t_{\mathrm{final}}/t_{\mathrm{start}})|$
$\triangleright$ number of iterations

13: for $k\leftarrow 1$ to $T$ do

14: $\delta_{x}:=-H(x)^{-1}\cdot\nabla f_{t}(x)$ $\triangleright$ $\delta_{x}$
is only used for analysis, not involve computation

15: $\triangleright$ $H(x)\in\mathbb{R}^{n\times n}$, $\nabla
f_{t}(x)\in\mathbb{R}^{n}$

16: $\widetilde{\delta}_{x}\leftarrow\textsc{StreamLS}(H(x),\nabla
f_{t}(x),\epsilon_{x}\cdot\epsilon_{\Phi})$ $\triangleright$ Algorithm 4

17: $\triangleright$ Compute any $\widetilde{\delta}_{x}$ that
$\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}\cdot\|\delta_{x}\|_{H(x)}$

18: $x^{\operatorname{new}}\leftarrow x+\widetilde{\delta}_{x}$

19: if $o=0$ then

20: $t^{\operatorname{new}}\leftarrow t\cdot(1-\epsilon_{t})$

21: else

22: $t^{\operatorname{new}}\leftarrow t\cdot(1+\epsilon_{t})$

23: end if

24: if ($o=0$ and $t<t_{\mathrm{final}}$) or ($o=1$ and
$t>t_{\mathrm{final}}$) then

25: break

26: end if

27: $t\leftarrow t^{\operatorname{new}}$, $x\leftarrow x^{\operatorname{new}}$

28: end for

29: return $x$

30:end procedure

##  7 Error analysis of IPM

The goal of this section is to present and prove Lemma 7.1.

Parameters | $\epsilon_{x}$ | $\epsilon_{t}$ | $\epsilon_{\Phi}$
---|---|---|---
Meaning | $\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}$ | $|t^{\operatorname{new}}/t-1|=\epsilon_{t}$ | $\Phi_{t}(x,x)\leq\epsilon_{\Phi}$
Value | $1/200$ | $1/(100\sqrt{m})$ | $1/100$
Table 2: Parameters occurred in Algorithm 3.

######  Lemma 7.1.

Given any matrix $A\in\mathbb{R}^{m\times n}$, $b\in\mathbb{R}^{m}$,
$c\in\mathbb{R}^{n}$. Let $x^{*}\in\mathbb{R}^{n}$ be the solution of the
following linear programming

$\min_{x\in\mathbb{R}^{n},Ax\geq b}c^{\top}x.$

If InteriorPointMethod is given
$A,b,c,t_{\mathrm{start}},x_{\mathrm{start}},t_{\mathrm{final}}$ that satisfy
$\Phi_{t_{\mathrm{start}}}(x_{\mathrm{start}},x_{\mathrm{start}})\leq\epsilon_{\Phi}$
(where $\epsilon_{\Phi}=1/100$, see Table 2), then it outputs an answer $x$
which is a _nearly-optimal solution_ :

$\displaystyle
c^{\top}x-c^{\top}x^{*}\leq\frac{m}{t_{\mathrm{final}}}\cdot(1+2\epsilon_{\Phi}).$

###### Proof.

Since our initial point $x_{\mathrm{start}}$ and $t_{\mathrm{start}}$ satisfy

$\Phi_{t_{\mathrm{start}}}(x_{\mathrm{start}},x_{\mathrm{start}})\leq\epsilon_{\Phi},$

we can apply Lemma 7.3 to get

$\Phi_{t_{\mathrm{final}}}(x,x)\leq\epsilon_{\Phi}.$

Applying Lemma 5.13 on $t_{\mathrm{final}}$ and $x$ with our choose of
$\epsilon_{\Phi}=1/100<1/10$, we get

$\displaystyle
c^{\top}x-c^{\top}x^{*}\leq\frac{m}{t_{\mathrm{final}}}\cdot(1+2\epsilon_{\Phi}).$

∎

The rest of this section is organized as follows: In Section 7.1, we state the
choices of our parameters. We will bound the potential function in Section
7.2. We bound the changes of the potential function when $t$ is moving in
Section 7.3. We bound the total movement of $x$ in Section 7.5. The detailed
proof of $x$ movement can be splitted into three parts: Section 7.6, Section
7.7 and Section 7.8.

###  7.1 Assumptions on parameters

######  Assumption 7.2.

Let $\epsilon_{x}\in(0,1/10)$, $\epsilon_{t}\in(0,1/\sqrt{m})$ and
$\epsilon_{\Phi}\in(0,1/10)$ denote three fixed parameters such that the
following inqueqality hold.

$\displaystyle
16(1+\epsilon_{t})\cdot\epsilon_{\Phi}^{2}+16\cdot\epsilon_{x}+\sqrt{m}\cdot\epsilon_{t}\leq\epsilon_{\Phi}.$

###  7.2 Bounding potential function $\Phi$

The goal of this section is to prove Lemma 7.3.

######  Lemma 7.3.

For each $t$, let
$\Phi_{t}:\mathbb{R}^{n}\times\mathbb{R}^{n}\rightarrow\mathbb{R}$ be defined
as Definition 5.7. Assume Assumption 7.2 holds. If the input
$x_{\mathrm{start}},t_{\mathrm{start}}$ satisfies

$\Phi_{t_{\mathrm{start}}}(x_{\mathrm{start}},x_{\mathrm{start}})\leq\epsilon_{\Phi},$

then for all iteration $k\in[T]$, we have

$\displaystyle\Phi_{t^{(k)}}(x^{(k)},x^{(k)})\leq\epsilon_{\Phi}.$

###### Proof.

We use proof by induction. In the basic case where $k=0$, we have
$t^{(0)}:=t_{\mathrm{start}}$ and $x^{(0)}:=x_{\mathrm{start}}$ so that the
condition holds by assumption.

When $k\geq 1$, for ease of notation, we define $x=x^{(k-1)}$, $t=t^{(k-1)}$,
$x^{\operatorname{new}}=x^{(k)}$, $t^{\operatorname{new}}=t^{(k)}$ and define
$\delta_{x}=H(x)^{-1}\nabla f_{t}(x)$ as in Definition 5.6.

First, from induction hypothesis, we have $\Phi_{t}(x,x)\leq\epsilon_{\Phi}$.
By definition of $\Phi$ (Definition 5.7), this means

$\displaystyle\|\delta_{x}\|_{H(x)}\leq\epsilon_{\Phi}.$

Next, we will show
$\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}$.

By calling $\textsc{StreamLS}(H(x),\nabla f_{t}(x),\epsilon_{x})$ (Line 16,
Algorithm 3) and applying Theorem 9.9, we have

$\displaystyle\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}\cdot\|\delta_{x}\|_{H(x)}.$

Define $b:=H(x)\delta_{x}$ and $v:=H(x)(\widetilde{\delta}_{x}-\delta_{x})$
and write $\widetilde{\delta}_{x}=H^{-1}\cdot(b+v)$. Note that we have
$\|v\|_{H(x)^{-1}}\leq\epsilon_{x}\cdot\|b\|_{H(x)^{-1}}$. Apply Lemma 5.14
and we get

$\widetilde{\delta}_{x}=\overline{H}^{-1}H(x)\delta_{x},$

where $\overline{H}\in\mathbb{R}^{n\times n}$ is a positive semi-definite
matrix that satisfies

$(1+20\epsilon_{x})^{-1}H(x)\preceq\overline{H}\preceq(1+20\epsilon_{x})H(x).$

Therefore, we have

$\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}=\|(\overline{H}^{-1}H-I)\delta_{x}\|_{H(x)}\leq
20\epsilon_{x}\|\delta_{x}\|_{H(x)}=20\epsilon_{x}\epsilon_{\Phi}\leq\epsilon_{x}.$

Finally, we have

$\displaystyle\Phi_{t^{\operatorname{new}}}(x^{\operatorname{new}},x^{\operatorname{new}})\leq$
$\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot\Phi_{t}(x^{\operatorname{new}},x^{\operatorname{new}})+|t^{\operatorname{new}}/t-1|\cdot\Psi(x^{\operatorname{new}},x^{\operatorname{new}})$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot\Phi_{t}(x^{\operatorname{new}},x^{\operatorname{new}})+|t^{\operatorname{new}}/t-1|\cdot\sqrt{m}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(1+\epsilon_{t})\cdot\Phi_{t}(x^{\operatorname{new}},x^{\operatorname{new}})+\sqrt{m}\cdot\epsilon_{t}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
16(1+\epsilon_{t})\cdot\Phi_{t}(x,x)^{2}+16\epsilon_{x}+\sqrt{m}\cdot\epsilon_{t}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
16(1+\epsilon_{t})\cdot\epsilon_{\Phi}^{2}+16\epsilon_{x}+\sqrt{m}\cdot\epsilon_{t}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\ \epsilon_{\Phi},$

where the first step is by Lemma 7.4, the second step is by Lemma 7.5, the
fourth step is by Lemma 7.6 since $\epsilon_{x}+\epsilon_{\Phi}\leq 1/4$, the
fifth step is by induction hypothesis, the last step is by Assumption 7.2. ∎

Lemma | Section | LHS | RHS
---|---|---|---
Lemma 7.3 | 7.2 | $\Phi_{t^{\operatorname{new}}}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})$ | $\epsilon_{\Phi}$
Lemma 7.4 | 7.3 | $\Phi_{t^{\operatorname{new}}}(x+\widetilde{\delta}_{x},x+\widetilde{\delta_{x}})$ | $(t^{\operatorname{new}}/t)\cdot\Phi_{t}(.)+(1-t^{\operatorname{new}}/t)\cdot\Psi(.)$
Lemma 7.5 | 7.4 | $\Psi(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})$ | $\sqrt{m}$
Lemma 7.6 | 7.5 | $\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})$ | $\Phi_{t}(x,x)^{2}$
Lemma 7.7 | 7.6 | $\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})$ | $\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})$
Lemma 7.8 | 7.7 | $\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})$ | $\Phi_{t}(x+\delta_{x},x+\delta_{x})$
Lemma 7.9 | 7.8 | $\Phi_{t}(x+\delta_{x},x+\delta_{x})$ | $\Phi_{t}(x,x)^{2}$
Table 3: Summary of movement of potential function.

###  7.3 Bounding the movement of $t$

The goal of this section is to prove Lemma 7.4.

######  Lemma 7.4.

Let $\Phi$ and $\Psi$ be defined in Definition 5.7 and Definition 5.8. Then,
for all positive $t,t^{\operatorname{new}}>0$ and feasible
$x\in\mathbb{R}^{n}$, we have

$\displaystyle\Phi_{t^{\operatorname{new}}}(x,x)\leq\frac{t^{\operatorname{new}}}{t}\cdot\Phi_{t}(x,x)+\left|t^{\operatorname{new}}/t-1\right|\cdot\Psi(x,x),$

###### Proof.

We have

$\displaystyle\nabla f_{t^{\operatorname{new}}}(x)=$
$\displaystyle\leavevmode\nobreak\ (t^{\operatorname{new}}\cdot c+g(x))$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot(t\cdot
c+g(x))+\left(1-\frac{t^{\operatorname{new}}}{t}\right)\cdot g(x)$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot\nabla
f_{t}(x)+\left(1-\frac{t^{\operatorname{new}}}{t}\right)\cdot g(x),$

where the first step follows from the definition of $f_{t}(x)$ (Definition
5.5), the second step follows from moving terms, and the last step follows
from the definition of $f_{t}(x)$ (Definition 5.5).

Finally, we can upper bound $\Phi_{t^{\operatorname{new}}}(x,x)$ as follows:

$\displaystyle\Phi_{t^{\operatorname{new}}}(x,x)=$
$\displaystyle\leavevmode\nobreak\ \|\nabla
f_{t^{\operatorname{new}}}(x)\|_{H(x)^{-1}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot\|\nabla
f_{t}(x)\|_{H(x)^{-1}}+\left(1-\frac{t^{\operatorname{new}}}{t}\right)\cdot\|g(x)\|_{H(x)^{-1}}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\frac{t^{\operatorname{new}}}{t}\cdot\Phi_{t}(x,x)+\left(1-\frac{t^{\operatorname{new}}}{t}\right)\cdot\Psi(x,x).$

Thus, we complete the proof. ∎

###  7.4 Upper bounding the potential function

The goal of this section is to prove Lemma 7.5.

######  Lemma 7.5.

Let function $\Psi:\mathbb{R}^{n}\times\mathbb{R}^{n}\rightarrow\mathbb{R}$ be
defined as Definition 5.8. For all feasible $x\in\mathbb{R}^{n}$, we have

$\displaystyle\Psi(x,x)\leq\sqrt{m}.$

###### Proof.

Let $z:=(H(x))^{-1}g(x)\in\mathbb{R}^{n}$ so that $\Psi(x,x)=\|z\|_{H(x)}$.

We can upper bound $\|z\|_{H(x)}^{2}$ as follows:

$\displaystyle\|z\|_{H(x)}^{2}=$ $\displaystyle\leavevmode\nobreak\
z^{\top}H(x)z$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
z^{\top}g(x)$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\sum_{i\in[m]}\frac{z^{\top}a_{i}}{s_{i}(x)}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\sqrt{m}\cdot\Big{(}\sum_{i\in[m]}\frac{(z^{\top}a_{i})^{2}}{s_{i}(x)^{2}}\Big{)}^{1/2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\sqrt{m}\cdot\Big{(}z^{\top}(\sum_{i\in[m]}\frac{a_{i}a_{i}^{\top}}{s_{i}(x)^{2}})z\Big{)}^{1/2}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\sqrt{m}\cdot(z^{\top}H(x)z)^{1/2}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ \sqrt{m}\cdot\|z\|_{H(x)},$

where the third step is by $g(x)=\sum_{i\in[m]}\frac{a_{i}}{s_{i}(x)}$
(Definition 5.4), the four step is by Cauchy-Schwarz, and the sixth step is by
$H(x)=\sum_{i\in[m]}\frac{a_{i}a_{i}^{\top}}{s_{i}(x)^{2}}$ (Definition 5.4).

Therefore, we obtain

$\displaystyle\Psi(x,x)=\|z\|_{H(x)}\leq\sqrt{m}.$

∎

###  7.5 Move both: final

The goal of this section is to prove Lemma 7.6

######  Lemma 7.6.

For any feasible $x\in\mathbb{R}^{n}$, let $\delta_{x}\in\mathbb{R}^{n}$ be
defined as Definition 5.6. Given
$\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}$ and
$\|\delta_{x}\|_{H(x)}\leq\epsilon_{1}$. If $\epsilon_{x}+\epsilon_{1}\leq
1/4$, then we have

$\displaystyle\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})\leq
16\cdot\Phi_{t}(x,x)^{2}+16\epsilon_{x}.$

###### Proof.

The proof is done by showing the following:

$\displaystyle\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})\leq$
$\displaystyle\leavevmode\nobreak\
\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})+16\epsilon_{x}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
16\cdot\Phi_{t}(x+\delta_{x},x+\delta_{x})+16\epsilon_{x}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ 16\cdot\Phi_{t}(x,x)^{2}+16\epsilon_{x},$

where the first step is by Lemma 7.7, the second step follows from Lemma 7.8,
and the third step follows from Lemma 7.9. ∎

###  7.6 Move both: part 1

The goal of this section is to prove Lemma 7.7.

######  Lemma 7.7.

Given $\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}$ and
$\|\delta_{x}\|_{H(x)}\leq\epsilon_{1}$. If $\epsilon_{x}+\epsilon_{1}\leq
1/4$, then we have

$\displaystyle\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})\leq\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})+16\epsilon_{x}.$

###### Proof.

We have

$\displaystyle\leavevmode\nobreak\
\Phi_{t}(x+\widetilde{\delta}_{x},x+\widetilde{\delta}_{x})-\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\ \|\nabla
f_{t}(x+\widetilde{\delta}_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}}-\|\nabla
f_{t}(x+\delta_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \|\nabla
f_{t}(x+\widetilde{\delta}_{x})-\nabla
f_{t}(x+\delta_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|g(x+\widetilde{\delta}_{x})-g(x+\delta_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}},$
(10)

where the first step is by the definition of $\Phi_{t}$ (Definition 5.7), and
the second step is by the triangle inequality.

Define
$\phi(t):=g\big{(}x+\delta_{x}+t\cdot(\widetilde{\delta}_{x}-\delta_{x})\big{)}$,
for $t\in[0,1]$. Then we have

$\displaystyle g(x+\widetilde{\delta}_{x})-g(x+\delta_{x})=\phi(1)-\phi(0).$

By mean value theorem, there exists $p\in[0,1]$ such that

$\displaystyle\phi(1)-\phi(0)=\phi^{\prime}(p)=H(x+\delta_{x}+p\cdot(\widetilde{\delta}_{x}-\delta_{x}))\cdot(\widetilde{\delta}_{x}-\delta_{x}).$

Since
$\|x+\delta_{x}+p\cdot(\widetilde{\delta}_{x}-\delta_{x})-x\|_{H(x)}\leq\|\delta_{x}\|_{H(x)}+p\cdot\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon+\epsilon_{1}\leq
1/4$, by Lemma 5.9 we have

$\displaystyle\frac{1}{4}H(x)\preceq
H(x+\delta_{x}+p\cdot(\widetilde{\delta}_{x}-\delta_{x}))\preceq 4H(x).$ (11)

Since
$\|x+\widetilde{\delta}_{x}-x\|_{H(x)}\leq\|\delta_{x}\|_{H(x)}+\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}+\epsilon_{1}\leq
1/4$, by Lemma 5.9 we have

$\displaystyle\frac{1}{4}H(x)\preceq H(x+\widetilde{\delta}_{x})\preceq
4H(x).$ (12)

Finally, we have

$\displaystyle\leavevmode\nobreak\
\|g(x+\widetilde{\delta}_{x})-g(x+\delta_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\big{\|}H(x+\delta_{x}+p\cdot(\widetilde{\delta}_{x}-\delta_{x}))\cdot(\widetilde{\delta}_{x}-\delta_{x})\big{\|}_{H(x+\widetilde{\delta}_{x})^{-1}}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
4\big{\|}H(x+\delta_{x}+p\cdot(\widetilde{\delta}_{x}-\delta_{x}))\cdot(\widetilde{\delta}_{x}-\delta_{x})\big{\|}_{H(x)^{-1}}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
16\big{\|}\widetilde{\delta}_{x}-\delta_{x}\big{\|}_{H(x)}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\ 16\epsilon_{x},$

where the first step is by Eq. (7.6), the second step is by Eq. (12), the
third step is by Eq. (11).

∎

###  7.7 Move both: part 2

The goal of this section is to prove Lemma 7.8.

######  Lemma 7.8.

Given $\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}$ and
$\|\delta_{x}\|_{H(x)}\leq\epsilon_{1}$. If $\epsilon_{x}+\epsilon_{1}\leq
1/4$, then we have

$\displaystyle\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})\leq
16\cdot\Phi_{t}(x+\delta_{x},x+\delta_{x}).$

###### Proof.

Since
$\|x+\widetilde{\delta}_{x}-x\|_{H(x)}\leq\|\delta_{x}\|_{H(x)}+\|\widetilde{\delta}_{x}-\delta_{x}\|_{H(x)}\leq\epsilon_{x}+\epsilon_{1}\leq
1/4$, by Lemma 5.9 we have

$\displaystyle\frac{1}{4}H(x)\preceq H(x+\widetilde{\delta}_{x})\preceq
4H(x).$

Since $\|x+\delta_{x}-x\|_{H(x)}\leq\epsilon_{1}\leq 1/4$, by Lemma 5.9 we
have

$\displaystyle\frac{1}{4}H(x)\preceq H(x+\delta_{x})\preceq 4H(x).$

Therefore,

$\displaystyle\Phi_{t}(x+\delta_{x},x+\widetilde{\delta}_{x})=$
$\displaystyle\leavevmode\nobreak\ \|\nabla
f_{t}(x+\delta_{x})\|_{H(x+\widetilde{\delta}_{x})^{-1}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ 4\|\nabla
f_{t}(x+\delta_{x})\|_{H(x)^{-1}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ 16\|\nabla
f_{t}(x+\delta_{x})\|_{H(x+\delta_{x})^{-1}}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\ 16\Phi_{t}(x+\delta_{x},x+\delta_{x}),$

where the first step is by definition of $\Phi$ (Definition 5.7), the second
step follows from $\frac{1}{4}H(x)\preceq H(x+\widetilde{\delta}_{x})\preceq
4H(x)$, the third step follows by $\frac{1}{4}H(x)\preceq
H(x+\delta_{x})\preceq 4H(x)$, the last step is by definition of $\Phi$
(Definition 5.7). ∎

###  7.8 Move both: part 3

The goal of this section is to prove Lemma 7.9.

######  Lemma 7.9.

$\displaystyle\Phi_{t}(x+\delta_{x},x+\delta_{x})\leq\Phi_{t}(x,x)^{2}$

###### Proof.

For simplicity, we define

$\displaystyle x^{\mathrm{tmp}}:=x+\delta_{x}\mathrm{\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ and\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ }z:=(H(x^{\mathrm{tmp}}))^{-1}\nabla
f_{t}(x^{\mathrm{tmp}}).$

First, we have

$\displaystyle\nabla f_{t}(x^{\mathrm{tmp}})$ $\displaystyle=t\cdot
c+g(x^{\mathrm{tmp}})$
$\displaystyle=-H(x)\cdot\delta_{x}-g(x)+g(x^{\mathrm{tmp}})$
$\displaystyle=-\sum_{i\in[m]}\frac{a_{i}a_{i}^{\top}}{s_{i}(x)^{2}}\delta_{x}+\sum_{i\in[m]}\frac{a_{i}}{s_{i}(x)}-\sum_{i\in[m]}\frac{a_{i}}{s_{i}(x^{\mathrm{tmp}})}$
$\displaystyle=\sum_{i\in[m]}\left(\frac{a_{i}a_{i}^{\top}}{s_{i}(x)s_{i}(x^{\mathrm{tmp}})}\cdot\delta_{x}-\frac{a_{i}a_{i}^{\top}}{s_{i}(x)^{2}}\cdot\delta_{x}\right)$
$\displaystyle=\sum_{i\in[m]}\frac{a_{i}\left(a_{i}^{\top}\delta_{x}\right)^{2}}{s_{i}(x)^{2}s_{i}(x^{\mathrm{tmp}})},$
(13)

where the first step is by definition of $f_{t}$ (Definition 5.5), the second
step is by $H(x)\cdot\delta_{x}=-t\cdot c-g(x)$, the third step is by
definition of $g(x)$ and $H(x)$ (Definition 5.4), and the fourth and last step
are by

$s_{i}(x^{\mathrm{tmp}})-s_{i}(x)=s_{i}(x+\delta_{x})-s_{i}(x)=a_{i}^{\top}(x+\delta_{x}-x)+(b-b)=a_{i}^{\top}\delta_{x}.$

As a result,

$\displaystyle\|z\|_{x^{\mathrm{tmp}}}^{2}=$
$\displaystyle\leavevmode\nobreak\ z^{\top}H(x^{\mathrm{tmp}})z$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\ z^{\top}\nabla
f_{t}(x^{\mathrm{tmp}})$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\sum_{i\in[m]}\frac{z^{\top}a_{i}\left(a_{i}^{\top}\delta_{x}\right)^{2}}{s_{i}(x)^{2}s_{i}(x^{\mathrm{tmp}})}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\Big{(}\sum_{i\in[m]}\frac{\left(z^{\top}a_{i}\right)^{2}}{s_{i}(x^{\mathrm{tmp}})^{2}}\Big{)}^{1/2}\cdot\sum_{i\in[m]}\frac{\left(a_{i}^{\top}\delta_{x}\right)^{2}}{s_{i}(x)^{2}}$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
(z^{\top}H(x^{\mathrm{tmp}})z)^{1/2}\cdot(\delta_{x}^{\top}H(x)\delta_{x})$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
\|z\|_{x^{\mathrm{tmp}}}\cdot\|\delta_{x}\|_{x}^{2},$

where the third step is by Eq. (7.8), the fourth step follows from Cauchy-
Schwarz inequality, the fifth step is by definition of $H(x)$ (Definition
5.4).

Therefore, we have that

$\displaystyle\Phi_{t}(x+\delta_{x},x+\delta_{x})=\|z\|_{x^{\mathrm{tmp}}}\leq\|\delta_{x}\|_{x}^{2}=\Phi_{t}(x,x)^{2}.$

∎

##  8 Pass complexity of IPM

The goal of this section is to prove Lemma 8.1. Recall that in the semi-
streaming model, the edge in the graph along with its weight is revealed one
by one in an adversarial order; the algorithm is permitted $\widetilde{O}(n)$
space and scans the edges for multiple passes; the edges order in different
passes need not to be the same.

######  Lemma 8.1 (Pass complexity of IPM).

Consider Algorithm 3 with input $A\in\mathbb{R}^{m\times n}$,
$b\in\mathbb{R}^{m}$, $c\in\mathbb{R}^{n}$,
$x_{\mathrm{start}}\in\mathbb{R}^{n}$, $t_{\mathrm{start}}\in\mathbb{R}$,
$t_{\mathrm{final}}\in\mathbb{R}$. Let $x\in\mathbb{R}^{n}$ be any feasible
point. Define $H(x)\in\mathbb{R}^{n\times n}$ and $s(x)\in\mathbb{R}^{m}$ as
the Hessian and slack variables (Definition 5.4, Definition 5.3). Suppose we
have the following conditions

  1. C.1

$\sum_{i\in[m]}a_{i}/s_{i}(x)$ can be computed in one pass in $f(n)$ space;

  2. C.2

the Hessian $H(x)\in\mathbb{R}^{n\times n}$ can be decomposed into

$H(x)=L_{G}(x)+D(x),$

where $L_{G}(x)\in\mathbb{R}^{n\times n}$ is the Laplacian of some graph $G$
with edge weight $w$ and $D(x)\in\mathbb{R}^{n\times n}$ is diagonal matrix;

  3. C.3

we can read the edge-weight pair $(e,w_{e})$ in one pass, and we can read the
diagonal of $D(x)\in\mathbb{R}^{n\times n}$ in one pass.

Then Algorithm 3 can be implemented in the streaming model with
$\widetilde{O}(n)+f(n)$ space and $\widetilde{O}(\sqrt{m})$ passes.

###### Proof.

Now consider the iterating part (Line 13-28). First we show $\nabla
f_{t}(x)\in\mathbb{R}^{n}$ can be computed in one pass.

By Definition 5.5,

$\nabla f_{t}(x)=c\cdot t-\sum_{i\in[m]}a_{i}/s_{i}(x),$

By condition C.1, $\nabla f_{t}(x)\in\mathbb{R}^{n}$ can be computed in one
pass in $f(n)$ space.

Since we have Conditions C.2 and C.3, we are able to first read $D$ in one
pass and store its entries in space $O(n)$, and then apply Theorem 9.9 to show
$\widetilde{\delta}_{x}$ can be computed in $\widetilde{O}(1)$ pass and
$\widetilde{O}(n)$ space.

All other computation can be done in $\widetilde{O}(n)$ space.

Since there are totally $T=O(\sqrt{m}\log(m/\epsilon_{\mathrm{ipm}}))$
iterations, the totally number of passes used is $\widetilde{O}(\sqrt{m})$.
Since we can reuse the space for $x$ and $t$, the total space used is
$\widetilde{O}(n)+f(n)$. ∎

##  9 SDD solver in the streaming model

In this section, we present an SDDM solver in the streaming model. Later in
Section A we reduce the problem of solving an SDD0 system to solving an SDDM
system and therefore give an SDD0 solver in the streaming model. The
definitions of these matrices can be found in Definition 3.3.

###  9.1 SDD and Laplacian systems

We need the definition of a spectral sparsifier, and a streaming algorithm for
computing it.

######  Definition 9.1 ($\delta$-spectral sparsifier).

Given a weighted undirected graph $G$ and a parameter $\delta>0$, an edge-
reweighted subgraph $H$ of $G$ is an _$\delta$ -spectral sparsifier_ of $G$
if101010We also say $L_{H}$ is an $\delta$-spectral sparsifier of $L_{G}$.

$(1-\delta)\cdot x^{\top}L_{G}x\leq x^{\top}L_{H}x\leq(1+\delta)\cdot
x^{\top}L_{G}x,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\
\forall x\in\mathbb{R}^{n},$

where $L_{G}$ and $L_{H}$ are the Laplacians of $G$ and $H$, respectively.

######  Lemma 9.2 ([KLM+17]).

Let $G$ be a weighted graph and $\delta\in(0,1)$ be a parameter. There exists
a streaming algorithm that takes $G$ as input, uses
$\delta^{-2}n\operatorname{poly}(\log n)$ space and $1$ pass, and outputs a
weighted graph $H$ with $\delta^{-2}n\operatorname{poly}(\log n)$ edges such
that with probability at least $1-1/\operatorname{poly}(n)$, $H$ is a
$\delta$-spectral sparsifier of $G$.

We will use the classic SDD solver in the sequential model, which is formally
described below.

######  Theorem 9.3 ([ST04]).

There is an algorithm which takes input an SDD0 matrix $A$, a vector
$b\in\mathbb{R}^{n}$, and a parameter $\epsilon\in(0,1/2)$, if there exists
$x^{*}\in\mathbb{R}^{n}$ such that $Ax^{*}=b$, then with probability
$1-1/\operatorname{poly}(n)$, the algorithm returns an $x\in\mathbb{R}^{n}$
such that $\|x-x^{*}\|_{A}\leq\epsilon\cdot\|x^{*}\|_{A}$ in

$\displaystyle\operatorname{nnz}(A)\cdot\operatorname{poly}(\log
n)\cdot\log(1/\epsilon)$

time. The returned $x$ is called an _$\epsilon$ -approximate solution_ to the
SDD0 system $Ax=b$.

###  9.2 The preconditioner

To prove that our SDDM solver (Algorithm 4) gives the desired accuracy, we
need the concept of a _preconditioner_ (and how to compute the preconditioner
of an SDDM matrix).

We define preconditioner as follows:

######  Definition 9.4 (Preconditioner).

For any positive definite matrix $A\in\mathbb{R}^{n\times n}$ and accuracy
parameter $\delta>0$, we say $P$ is a $\delta$-precondioner if

$\displaystyle\|P^{-1}Ax-x\|_{A}\leq\delta\cdot\|x\|_{A},\leavevmode\nobreak\
\leavevmode\nobreak\ \leavevmode\nobreak\ \forall x\in\mathbb{R}^{n}.$

######  Lemma 9.5.

Let $A$ be an SDDM matrix (Definition 3.3) and let $A=L_{G}+D$ where $L_{G}$
is a Laplacian matrix of graph $G$ and $D\succ 0$ is a diagonal matrix with
non-negative entries.111111By Fact 3.4, such decomposition always exists. For
any $\delta\in(0,1/2)$, if $H$ is a $\delta$-spectral sparsifier of $G$, and
we define $P:=L_{H}+D$. Then $P$ satisfies the following two conditions:

  * •

$P$ is a $(2\delta)$-preconditioner(Definition 9.4) of $A$,

  * •

$P\succ 0$.

###### Proof.

The lemma clearly holds for $x=\mathbf{0}_{n}$, so assume not in the
following. By Definition 9.1,

$(1+\delta)x^{\top}L_{G}x\geq x^{\top}L_{H}x\geq(1-\delta)x^{\top}L_{G}x\geq
0,$

then $L_{H}\succeq 0$ as it must be symmetric. Since $D\succ 0$, we have that
$A\succ 0$, $P\succ 0$, and

$(1+\delta)x^{\top}(L_{G}+D)x\geq
x^{\top}(L_{H}+D)x\geq(1-\delta)x^{\top}(L_{G}+D)x>0,$

which is $(1+\delta)A\succeq P\succeq(1-\delta)A$. So we obtain

$\frac{1}{1-\delta}A^{-1}\succeq P^{-1}\succeq\frac{1}{1+\delta}A^{-1},$

which, by $\delta\in(0,1/2)$, implies

$(1+2\delta)A^{-1}\succeq P^{-1}\succeq(1-2\delta)A^{-1}.$ (14)

Since all matrices in Eq. (14) are positive definite, we can left multiply
$A^{1/2}$ and right multiply $A^{1/2}$ on Eq. (14) to get

$(1+2\delta)I\succeq A^{1/2}P^{-1}A^{1/2}\succeq(1-2\delta)I.$ (15)

Subtracting $I$ from Eq. (15) then left multiplying $x^{\top}$ and right
multiplying $x$, we obtain

$\displaystyle 2\delta x^{\top}x\geq
x^{\top}\underbrace{({A}^{1/2}P^{-1}{A}^{1/2}-I)}_{:=M}x\geq-2\delta
x^{\top}x.$ (16)

Since Eq. (16) holds for any $x\in\mathbb{R}^{n}$ and $M$ is symmetric, using
spectral theorem, we get

$2\delta\geq\lambda_{n}(M)\geq\lambda_{1}(M)\geq-2\delta.$

Next, we have that

$\displaystyle\|Mx\|_{2}=\max\\{|\lambda_{n}(M)|,|\lambda_{1}(M)|\\}\leq
2\delta\|x\|_{2}.$

Let $y:={A}^{-1/2}x$. As $\text{rank}({A}^{-1/2})=n$, we get that for any
$y\in\mathbb{R}^{n}$,

$2\delta\|{A}^{1/2}y\|_{2}\geq\|M{A}^{1/2}y\|_{2}=\|{A}^{1/2}P^{-1}Ay-{A}^{1/2}y\|_{2}.$
(17)

Finally, rewriting both sides of Eq. (17) by the definition of matrix norm, we
obtain

$2\delta\|y\|_{A}\geq\|P^{-1}Ay-y\|_{A}$

for any $y\in\mathbb{R}^{n}$, giving the lemma. ∎

###  9.3 An iterative solver

In this section, we present a streaming SDDM solver (Algorithm 4) that takes
$\widetilde{O}(n)$ space and $O(\log(1/\epsilon)/\log\log n)$ passes for
approximately solving SDDM system $Ax=b$ with error parameter
$\epsilon\in(0,1/10)$.

Algorithm 4 $\textsc{StreamLS}(A=L_{G}+D,b,\epsilon)$. Return an
$\epsilon$-approx solution to $A^{-1}b$

1:procedure StreamLS($A=L_{G}+D\in\mathbb{R}^{n\times n}$,
$b\in\mathbb{R}^{n}$, $\epsilon\in\mathbb{R}$)

2: $\triangleright$ Note that $L_{G}$ is from input stream. $D$ is diagonal
matrix stored.

3: $\delta\leftarrow 1/\log n$

4: $\widetilde{\delta}\leftarrow\delta/2$

5: $T\leftarrow O(\max\\{1,(\log(1/\epsilon))/{\log\log n}\\})$
$\triangleright$ Number of iterations

6: Compute and store a $\delta$-spectral sparsifier $H$ of $L_{G}$
$\triangleright$ Use one pass, Lemma 9.2

7: Compute a $\delta$-preconditioner of $A$ as $P:=L_{H}+D$

8: $r_{0}\leftarrow b$, $x_{0}\leftarrow\mathbf{0}_{n}$ $\triangleright$
$r_{0},x_{0}\in\mathbb{R}^{n}$

9: for $t\leftarrow 0$ to $T-1$ do

10: Compute a $\widetilde{\delta}$-approximate solution $y_{t}$ to $Py=r_{t}$
by an SDD0 solver and store $y_{t}$ in the memory $\triangleright$ Theorem 9.3

11: $r_{t+1}\leftarrow r_{t}-Ay_{t}$

12:

13: $x_{t+1}\leftarrow x_{t}+y_{t}$

14: end for

15: return $x_{T}$.

16:end procedure

###  9.4 An iterative solver: the space

The goal of this section is to prove Lemma 9.6.

######  Lemma 9.6.

Let $A=L_{G}+D\in\mathbb{R}^{n\times n}$ be an SDDM matrix where
$L_{G}\in\mathbb{R}^{n\times n}$ is the Laplacian matrix of graph $G$ with
weight $w$, and $D\in\mathbb{R}^{n\times n}$ is a diagonal matrix. If we can
read all edge-weight pairs $(e,w_{e})$ in one pass, and if we can read the
diagonal of matrix $D$ in one pass, then $\textsc{StreamLS}(A,b,\epsilon)$
takes $O(\max\\{1,\log(1/\epsilon)/\log\log n\\})$ passes and
$\widetilde{O}(n)$ space.121212The algorithm can be implemented in the
standard RAM model with finite precision by introducing an $\widetilde{O}(1)$
factor in the encoding, which translates to a multiplicative factor of
$\widetilde{O}(1)$ in the space [ST14].

###### Proof.

By Lemma 9.2 and $\delta=1/\log n$, the $\delta$-spectral sparsifier $H$ has
$\widetilde{O}(n)$ edges and can be computed in $1$ pass and
$\widetilde{O}(n)$ space with probability $1-1/\operatorname{poly}(n)$.
Therefore, computing the $\delta$-preconditioner $P$ also takes $1$ pass and
$\widetilde{O}(n)$ space. Note that $P\succ 0$ and thus the system $Py=r_{t}$
always has a solution.

It remains to prove that each iteration of Lines 9-14 takes $1$ pass and
$\widetilde{O}(n)$ space. Since any iteration $t$ only needs the vectors
subscripted by $t$ and $t+1$, we can reuse the space such that the total space
is $\widetilde{O}(n)$.

Since $\text{nnz}(P)=\widetilde{O}(n)$ and

$\widetilde{\delta}=\delta/2=1/(2\log n),$

in Line 10, by Theorem 9.3, with probability $1-1/\operatorname{poly}(n)$ a
$\widetilde{\delta}$-approximate solution $y_{t}$ can be found in
$\widetilde{O}(n\log(1/\widetilde{\delta}))=\widetilde{O}(n)$ time, and
therefore in $\widetilde{O}(n)$ space. Note that this step does not read the
stream.

In Line 11, computing $r_{t+1}$ requires computing $Ay_{t}$, which is done by
reading the stream of $L_{G}$ and $D$ for $1$ pass and multiplying the
corresponding entries and adding up to the corresponding coordinate. All
vectors are in $\mathbb{R}^{n}$, so the total space used in each iteration is
$\widetilde{O}(n)$. The lemma follows immediately from

$\displaystyle T=O(\max\\{1,\log(1/\epsilon)/\log\log n\\})$

and a union bound. ∎

###  9.5 An iterative solver: the accuracy

The goal of this section is to prove Lemma 9.7 (which shows that the solver
iteratively improves the accuracy of the solution).

######  Lemma 9.7.

Let $\delta\in(0,1/10)$. For any $t\in[0,T]$,
$\|x_{t}-A^{-1}b\|_{A}\leq(4\delta)^{t}\cdot\|A^{-1}b\|_{A}$.

###### Proof.

The proof is by an induction on $t$. In the basic case of $t=0$, we have
$x_{t}=0$ so the statement clearly holds. Assuming the lemma holds for $t$, we
prove the inductive step for $t+1$.

Since $y_{t}$ is a $\widetilde{\delta}$-approximate solution (Line 10,
Algorithm 4), by Theorem 9.3 we have

$\|y_{t}-P^{-1}r_{t}\|_{P}\leq\widetilde{\delta}\cdot\|P^{-1}r_{t}\|_{P}.$

By definition of the matrix norm, this becomes

$(y_{t}-P^{-1}r_{t})^{\top}P(y_{t}-P^{-1}r_{t})\leq\widetilde{\delta}^{2}\cdot(P^{-1}r_{t})^{\top}P(P^{-1}r_{t}).$
(18)

Since $P$ is a $\delta$-preconditioner, assuming $y_{t}\neq P^{-1}r_{t}$ and
applying Definition 9.1 on both sides of Eq. (18), we have

$(1-\delta)\cdot(y_{t}-P^{-1}r_{t})^{\top}A(y_{t}-P^{-1}r_{t})\leq(1+\delta)\widetilde{\delta}^{2}\cdot(P^{-1}r_{t})^{\top}A(P^{-1}r_{t}).$
(19)

If $y_{t}=P^{-1}r_{t}$, then Eq. (19) also holds since the right-hand side of
Eq. (19) is non-negative due to $A\succeq 0$. Note that Eq. (19) implies

$\displaystyle\|y_{t}-P^{-1}r_{t}\|_{A}\leq$
$\displaystyle\leavevmode\nobreak\
\widetilde{\delta}\sqrt{(1+\delta)/(1-\delta)}\cdot\|P^{-1}r_{t}\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
2\widetilde{\delta}\cdot\|P^{-1}r_{t}\|_{A}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \delta\cdot\|P^{-1}r_{t}\|_{A},$ (20)

where the second step follows from $(1+\delta)/(1-\delta)\leq 4$ when
$\delta\in(0,1/2)$, and the last step follows from
$\widetilde{\delta}=\delta/2$.

Before continuing, we observe the following, which easily follows from the
update rule and an induction on $t$:

$r_{t}=r_{0}-A\sum_{i=0}^{t-1}y_{i}=r_{0}-A(x_{t}-x_{0})=b-Ax_{t}.$ (21)

Using Eq. (21), we bound the left-hand side of Eq. (9.5) from above by the
following:

$\displaystyle\leavevmode\nobreak\ \|y_{t}-P^{-1}r_{t}\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\delta\cdot\|P^{-1}(b-Ax_{t})\|_{A}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\delta\cdot\|P^{-1}A(A^{-1}b-x_{t})-(A^{-1}b-x_{t})+(A^{-1}b-x_{t})\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\delta\cdot\left(\|P^{-1}A(A^{-1}b-x_{t})-(A^{-1}b-x_{t})\|_{A}+\|A^{-1}b-x_{t}\|_{A}\right)$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\delta\cdot\left(2\delta\|A^{-1}b-x_{t}\|_{A}+\|A^{-1}b-x_{t}\|_{A}\right)$
$\displaystyle=$ $\displaystyle\leavevmode\nobreak\
2\delta\cdot\|A^{-1}b-x_{t}\|_{A},$ (22)

where the third step follows from the triangle inequality, the fourth step
follows from Lemma 9.5, the last step follows from $\delta\in(0,1/2)$.

Finally, we are ready to prove the inductive step:

$\displaystyle\|x_{t+1}-A^{-1}b\|_{A}=$ $\displaystyle\leavevmode\nobreak\
\|x_{t}-A^{-1}b+y_{t}\|_{A}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|(x_{t}-A^{-1}b+P^{-1}r_{t})+(y_{t}-P^{-1}r_{t})\|_{A}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
\|x_{t}-A^{-1}b+P^{-1}r_{t}\|_{A}+\|y_{t}-P^{-1}r_{t}\|_{A}$ $\displaystyle=$
$\displaystyle\leavevmode\nobreak\
\|x_{t}-A^{-1}b+P^{-1}b-P^{-1}Ax_{t}\|_{A}+\|y_{t}-P^{-1}r_{t}\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\|P^{-1}A(A^{-1}b-x_{t})-(A^{-1}b-x_{t})\|_{A}+2\delta\cdot\|A^{-1}b-x_{t}\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
2\delta\cdot\|A^{-1}b-x_{t}\|_{A}+2\delta\cdot\|{A}^{-1}b-x_{t}\|_{A}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(4\delta)^{t+1}\cdot\|A^{-1}b\|_{A},$

where the first step follows from $x_{t+1}=x_{t}+y_{t}$, the third step
follows from the triangle inequality, the fourth step follows from Eq. (21),
the fifth step follows from Eq. (9.5), the sixth step follows from Lemma 9.5,
and the last step follows from the induction hypothesis, completing the proof.
∎

###  9.6 Main result

The goal of this section is to prove Theorem 9.9

######  Lemma 9.8.

Let $\epsilon\in(0,1)$ and $b\in\mathbb{R}^{n}$. Let
$A=L_{G}+D\in\mathbb{R}^{n\times n}$ be an SDDM matrix where
$L_{G}\in\mathbb{R}^{n\times n}$ is the Laplacian matrix of graph $G$ with
weight $w$, and $D\in\mathbb{R}^{n\times n}$ is a diagonal matrix. If we can
read all edge-weight pairs $(e,w_{e})$ in one pass, and if we can read the
diagonal of matrix $D$ in one pass, then with probability
$1-1/\operatorname{poly}(n)$, $\textsc{StreamLS}(A,b,\epsilon)$ returns an
$\epsilon$-approximate solution $x$, i.e.,

$\|x-A^{-1}b\|_{A}\leq\epsilon\cdot\|A^{-1}b\|_{A}.$

in $O(\max\\{1,\log(1/\epsilon)/\log\log n\\})$ passes and $\widetilde{O}(n)$
space.

###### Proof.

It follows by Lemma 9.6, our choices of $\delta,T$, and Lemma 9.7. ∎

By Lemma 9.8 and the reduction in Section A, we obtain our main result.

######  Theorem 9.9.

There is a streaming algorithm which takes input an SDD0 matrix $A$, a vector
$b\in\mathbb{R}^{n}$, and a parameter $\epsilon\in(0,1)$, if there exists
$x^{*}\in\mathbb{R}^{n}$ such that $Ax^{*}=b$, then with probability
$1-1/\operatorname{poly}(n)$, the algorithm returns an $x\in\mathbb{R}^{n}$
such that $\|x-x^{*}\|_{A}\leq\epsilon\cdot\|x^{*}\|_{A}$ in
$O(\max\\{1,\log(1/\epsilon)/\log\log n\\})$ passes and $\widetilde{O}(n)$
space.

##  10 Minimum vertex cover

In this section, our goal is that given an linear programming form of
fractional minimum vertex cover, output the set of tight constraints of some
optimal solution. This set is crucial for us to turn the optimal dual solution
into an optimal primal solution via complementary slackness (Theorem 11.6).

The basic idea is to run interior point method (procedure InteriorPointMethod,
Section 6, Algorithm 3) to get a near-optimal solution. However, the hard part
is that even if we get a near-optimal solution, we are still far from figuring
out all tight constraints of “some” optimal solution, e.g., when $\langle
x,c\rangle$ is close to $\langle x^{*},c\rangle$, $x$ is not necessarily close
to $x^{*}$. The key observation is that our LP is integral and does not have
large bit complexity, which means once we get a near-optimal solution, we can
align to some optimal solution. Technically, here we use the isolation lemma
for the second time (Lemma 10.2).

The goal of this section is to prove Lemma 10.1 (Correctness) and Lemma 10.4
(Pass complexity). As a by-product, we also give a formal theorem on solving
fractional minimum vertex cover.

Statement | Section | Where we use | Comments
---|---|---|---
Lemma 10.1 | 10.2 | Theorem 11.3 | Main lemma, correctness
Lemma 10.4 | 10.3 | Theorem 11.2 | Main lemma, running time
Lemma 10.2 | 10.2 | Lemma 10.1 | From nearly-optimal solution to exact solution
Lemma 10.3 | 10.2 | Lemma 10.1 | Bound the eigenvalue of Hessian
Lemma 10.5 | 10.3 | Lemma 10.4 | Auxiliary lemma for running time
Lemma 10.6 | 10.4 | Lemma 10.3, 10.4 | Prove Hessian matrix is SDDM matrix
Table 4: The structure of Section 10.

###  10.1 Our algorithm

Algorithm 5 Minimum vertex cover

1:procedure
MinimumVertexCover($G=(V_{L},V_{R},E),b_{1}\in\operatorname*{\mathbb{Z}}^{m},c_{1}\in\operatorname*{\mathbb{Z}}^{n},n,m$)
$\triangleright$ Lemma 10.1

2: Let $A_{1}$ be the signed edge-vertex incident matrix of $G$ with direction
$V_{L}$ to $V_{R}$

3: $L\leftarrow$ the bit complexity of $A_{1},b_{1},c_{1}$

4: Modify $A_{1},b_{1},c_{1}$ by adding constraints $x_{V_{L}}\geq 0$ and
$x_{V_{R}}\leq 0$ to get $A_{2},b_{2},c_{2}$

5: Modify $A_{2},b_{2},c_{2}$ according to Lemma 10.2 and get
$A_{3},b_{3},c_{3}$ $\triangleright$ $A_{2}=A_{3}$, $b_{2}=b_{3}$

6: $x_{\mathrm{init}}\leftarrow
2^{L}\cdot(\mathbf{1}_{V_{L}}-\mathbf{1}_{V_{R}})$

7: $t_{\mathrm{init}}\leftarrow 1$

8: $c_{\mathrm{init}}\leftarrow-\nabla g(x_{\mathrm{init}})$

9: $\epsilon_{\Phi}\leftarrow 1/100$

10: $t_{1}\leftarrow\epsilon_{\Phi}\cdot(m2^{4L+10})^{-1}$

11:
$x_{\mathrm{tmp}}\leftarrow\textsc{InteriorPointMethod}(A_{3},b_{3},c_{\mathrm{init}},t_{\mathrm{init}},x_{\mathrm{init}},t_{1},0,\epsilon_{\Phi}/4)$
$\triangleright$ Algorithm 3

12: $t_{2}\leftarrow nm2^{3L+10}$

13:
$x\leftarrow\textsc{InteriorPointMethod}(A_{3},b_{3},c_{3},t_{1},x_{\mathrm{tmp}},t_{2},\epsilon_{\Phi},1)$
$\triangleright$ Algorithm 3

14: Use Lemma 10.2 to turn $x$ into a set $S$ of tight constraints on LP of
$(A_{1},b_{1},c_{1})$

15: $\triangleright$ $S\subseteq[m+2n]$

16: return $S\cap[m]$

17:end procedure

Description. The constraints of minimum vertex cover is $Ax\geq
b,\leavevmode\nobreak\ x\geq 0$ while our IPM algorithm only accepts the form
$Ax\geq b$. In order to call IPM, at the beginning we transform
$A_{1},b_{1},c_{1}$ to $A_{2},b_{2},c_{2}$. Then by Lemma 10.2, we do slight
perturbation over $A_{2},b_{2},c_{2}$ to get $A_{3},b_{3},c_{3}$. Then the
algorithm executes the standard two IPM walks. The first walk starts from a
far point on the central path which is ensured to be feasible and then goes
along the central path to the analytic center. Then we switch the central path
which correspond to the true LP we want to solve, starting from this analytic
center and solve it. Basically, the first walk is to find a good initial point
for the second walk.

###  10.2 Correctness of Algorithm 5

The goal of this section is to prove Lemma 10.1.

######  Lemma 10.1 (Correctness of Algorithm 5).

In Algorithm 5, let $A=A_{1}$ be signed edge-vertex incident matrix of $G$
with direction $V_{L}$ to $V_{R}$. Let $[n]=V_{L}\cup V_{R}$. Let $b=b_{1}$
and $c=c_{1}$ be the input. If the linear programming

$\displaystyle\min_{x\in\mathbb{R}^{n}}$ $\displaystyle\leavevmode\nobreak\
c^{\top}x$ (23) $\displaystyle\mathrm{s.t.}$
$\displaystyle\leavevmode\nobreak\ Ax\geq b$
$\displaystyle\leavevmode\nobreak\ x_{v}\geq 0,\forall v\in V_{L}$
$\displaystyle\leavevmode\nobreak\ x_{v}\leq 0,\forall v\in V_{R}$

is both feasible and bounded, then with probability at least $1/2$, Algorithm
5 returns a set of tight constraints on some optimal solution to the LP
Eq.(23).

###### Proof.

First we need to show these two calls into InteriorPointMethod satisfy the
initial condition that $x_{\mathrm{start}},t_{\mathrm{start}}$ is a good start
point.

The first call In the first call (Line 11), we have

$\displaystyle\|t_{\mathrm{init}}\cdot c_{\mathrm{init}}+\nabla
g(x_{\mathrm{init}})\|_{H(x_{\mathrm{init}})^{-1}}=\|0\|_{H(x_{\mathrm{init}})^{-1}}=0\leq\epsilon_{\Phi}/4,$

where the first step is by definition of $c_{\mathrm{init}}$.

Since the algorithm InteriorPointMethod ends up in parameter
$x_{\mathrm{tmp}}$ and $t_{1}$, by Lemma 7.3, we have

$\displaystyle\|t_{1}\cdot c_{\mathrm{init}}+\nabla
g(x_{\mathrm{tmp}})\|_{H(x_{\mathrm{tmp}})^{-1}}\leq\epsilon_{\Phi}/4.$ (24)

Note that $t_{1}$ and $x_{\mathrm{tmp}}$ is the input to the second call.

The second call In the second call (Line 13), the desired term can be upper
bounded by

$\displaystyle\leavevmode\nobreak\ \|t_{1}\cdot c_{3}+\nabla
g(x_{\mathrm{tmp}})\|_{H(x_{\mathrm{tmp}})^{-1}}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ \|t_{1}\cdot c_{\mathrm{init}}+\nabla
g(x_{\mathrm{tmp}})\|_{H(x_{\mathrm{tmp}})^{-1}}+\|t_{1}\cdot
c_{\mathrm{init}}-t_{1}\cdot c_{3}\|_{H(x_{\mathrm{tmp}})^{-1}}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\epsilon_{\Phi}/4+t_{1}\cdot\|c_{\mathrm{init}}-c_{3}\|_{H(x_{\mathrm{tmp}})^{-1}}$

where the first step is by triangle inequality, the second step is from Eq.
(24).

Now let’s bound the term
$\|c_{\mathrm{init}}-c_{3}\|_{H(x_{\mathrm{tmp}})^{-1}}$.

$\displaystyle\|c_{\mathrm{init}}-c_{3}\|_{H(x_{\mathrm{tmp}})^{-1}}\leq$
$\displaystyle\leavevmode\nobreak\
\|c_{\mathrm{init}}-c_{3}\|_{2}\cdot\lambda_{\min}(H(x_{\mathrm{tmp}}))^{-1/2}$
$\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
\|c_{\mathrm{init}}-c_{3}\|_{2}\cdot 2^{L+3}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\
(\|c_{\mathrm{init}}\|_{2}+\|c_{3}\|_{2})\cdot 2^{L+3}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ (\sqrt{m}2^{L+3}+\|c_{3}\|_{2})\cdot
2^{L+3}$ $\displaystyle\leq$ $\displaystyle\leavevmode\nobreak\
(\sqrt{m}2^{L+3}+m2^{3L+4})\cdot 2^{L+3}$ $\displaystyle\leq$
$\displaystyle\leavevmode\nobreak\ m2^{4L+8},$

where the first step is by
$\|x\|_{H}=\sqrt{x^{\top}Hx}\leq\|x\|_{2}\cdot\sqrt{\lambda_{\max}(H)}$ and
$\lambda_{\max}(H)=\lambda_{\min}(H^{-1})^{-1}$, the second step is by Lemma
10.3, the fourth step is by $c_{\mathrm{init}}=-\nabla g(x_{\mathrm{init}})$,
$x_{\mathrm{init}}=2^{L}(\mathbf{1}_{V_{L}}-\mathbf{1}_{V_{R}})$ and the
definition of $g$ (Definition 5.4), the fifth step is by definition of $c_{3}$
(in Lemma 10.2).

By our choice of $t_{1}:=\epsilon_{\Phi}\cdot(m2^{4L+10})^{-1}$ (Line 10), we
finally have

$\|t_{1}\cdot c_{3}+\nabla
g(x_{\mathrm{tmp}})\|_{H(x_{\mathrm{tmp}})^{-1}}\leq\epsilon_{\Phi}/4+\epsilon_{\Phi}/4\leq\epsilon_{\Phi}.$

Let $\mathrm{OPT}$ be the optimal solution to the linear programming

$\displaystyle\min_{x\in\mathbb{R}^{n},A_{3}x\geq b_{3}}c_{3}^{\top}x.$

By $t_{2}:=nm2^{3L+10}$ and Lemma 7.1, we know our solution $x$ (Line 13) is
feasible and has value
$c_{3}^{\top}x-\mathrm{OPT}\leq\frac{m}{t_{2}}\cdot(1+2\epsilon_{\Phi})\leq
n^{-1}2^{-3L-2}$. Now, we can apply Lemma 10.2 to show that after one matrix
vector multiplication, with probability at least $1/2$, we can output the
tight constraints $S$ of a basic feasible optimal solution of

$\displaystyle\min_{x\in\mathbb{R}^{n},A_{2}x\geq b_{2}}c_{2}^{\top}x.$ (25)

Note that this LP is exactly the same as Eq. (23) by our construction on
$A_{2},b_{2},c_{2}$ (Line 4).

∎

######  Lemma 10.2 (Lemma 43 of [LS14]).

Given a feasible and bounded linear programming

$\displaystyle\leavevmode\nobreak\ \min_{x\in\mathbb{R}^{n}}c^{\top}x$ (26)
$\displaystyle\mathrm{s.t.\leavevmode\nobreak\ }$
$\displaystyle\leavevmode\nobreak\ Ax\geq b,$

where $A\in\operatorname*{\mathbb{Z}}^{m\times n}$,
$b\in\operatorname*{\mathbb{Z}}^{m}$, $c\in\operatorname*{\mathbb{Z}}^{n}$ all
having integer coefficient. Let $L$ be the bit complexity of Eq. (26).

Let $r\in\operatorname*{\mathbb{Z}}^{n}$ be chosen uniformly at random from
the integers $\\{-2^{L+1}n,\cdots,2^{L+1}n\\}$. Consider the perturbed linear
programming

$\displaystyle\leavevmode\nobreak\ \min_{x\in\mathbb{R}^{n}}(2^{2L+3}n\cdot
c+r)^{\top}x$ (27) $\displaystyle\mathrm{s.t.\leavevmode\nobreak\ }$
$\displaystyle\leavevmode\nobreak\ Ax\geq b.$

Then with probability at least $1/2$ over the randomness on
$r\in\operatorname*{\mathbb{Z}}^{n}$, we have the following.

Let $\mathrm{OPT}$ be defined as the optimal value of linear programming Eq.
(27). Let $x$ be any feasible solution for Eq. (27) with objective value less
than $\mathrm{OPT}+n^{-1}2^{-3L-2}$, then we can find the tight constraints of
a basic feasible optimal solution of Eq. (26) using one matrix vector
multiplication with $A$. Moreover, we have $\|x-x^{*}\|_{\infty}\leq 1/n$,
where $x^{*}$ is the unique optimal solution for Eq. (27). Plus, the bit
complexity of Eq. (27) is at most $3L+\log(8n)$.

######  Lemma 10.3.

Let $H(x)$ be defined as in Lemma 10.6, then we have

$\displaystyle\lambda_{\min}(H(x))\geq 2^{-2L-6}.$

###### Proof.

We can lower bound $\lambda_{\min}(H(x))$ in the following sense,

$\displaystyle\lambda_{\min}(H(x))\geq$ $\displaystyle\leavevmode\nobreak\
\lambda_{\min}(D(x))$ $\displaystyle\geq$ $\displaystyle\leavevmode\nobreak\
\min_{i\in[n]}D(x)_{i,i}$ $\displaystyle\geq$
$\displaystyle\leavevmode\nobreak\ 2^{-2L-6}$

where the first step is by Lemma 10.6 that $H(x)=L(x)+D(x)$ and both $L(x)$
and $D(x)$ are positive semi-definite matrices, the second step is by the fact
that $D(x)$ is diagonal matrix, the third step is by Lemma 10.6 that
$D(x)_{i,i}=s_{i+m}(x)^{-2}+s_{i+m+n}(x)^{-2}$ and $s(x):=Ax-b$ is bounded by
$\|s\|_{\infty}\leq 2^{L+3}$. ∎

###  10.3 Pass complexity of Algorithm 5

######  Lemma 10.4 (Pass complexity of Algorithm 5).

Suppose there’s an oracle $\mathcal{I}$ running in $f(n)$ space that can
output $(b_{1})_{i}$ given any $i\in[m]$. If the input satisfies
$\|b\|_{\infty},\|c\|_{\infty}$ is polynomially bounded, then Algorithm 5 can
be implemented in the streaming model within $\widetilde{O}(n)+f(n)+|S|$ space
and $\widetilde{O}(\sqrt{m})$ passes.

###### Proof.

In Line 2, we defined $A_{1}$ but never explicitly compute and store $A_{1}$
in memory.

In Line 3, we let $L\leftarrow O(\log n)$ as the upper bound on the bit
complexity. This cost one pass. Indeed, Lemma 3.5 implies that
$|\mathrm{d}_{\mathrm{max}}(A_{1})|\leq 1$. So $L$ is the upper bound on the
bit complexity.

In Line 4 and Line 5, we defined $A_{2},b_{2},c_{2},A_{3},b_{3},c_{3}$. This
doesn’t involve computation.

In Lines 6, 7, 8, 9, 10, 12, they are done locally in memory.

In Line 11 and Line 13, we will show the following three things in order to
call Lemma 8.1.

  1. 1.

$\sum_{i\in[m]}a_{i}/s_{i}(x)$ can be computed in one pass in $h(n)$ space;

  2. 2.

the Hessian $H(x)\in\mathbb{R}^{n\times n}$ can be decomposed into

$H(x)=L_{G}(x)+D(x),$

where $L_{G}(x)\in\mathbb{R}^{n\times n}$ is the Laplacian of some graph $G$
with edge weight $w$ and $D(x)\in\mathbb{R}^{n\times n}$ is diagonal matrix;

  3. 3.

we can read the edge-weight pair $(e,w_{e})$ in one pass, and we can read the
diagonal of $D(x)$ in one pass.

Part 1: By Lemma 10.5 and our assumption of the existence of such oracle
$\mathcal{I}$, given a feasible $x\in\mathbb{R}^{n}$, we can output a stream
of all edges $e_{i}=(u,v)$ together with its slack
$s_{i}(x)=a_{i}^{\top}x-b_{i}$ in one pass within $O(n)+f(n)$ space. So
$\sum_{i\in[m]}a_{i}/s_{i}(x)$ can be computed in one pass in $h(n)=O(n)+f(n)$
space.

Part 2: This follows by Lemma 10.6.

Part 3: By Lemma 10.6, the edge weight is $w_{i}=s_{i}(x)^{-2}$, and
$D(x)_{i,i}=s_{i+m}(x)^{-2}+s_{i+m+n}^{-2}$. So Part 3 is satisfied by Lemma
10.5.

By Lemma 8.1, Line 11 and Line 13 can be implemented in streaming model within
$\widetilde{O}(n)+f(n)$ space and $\widetilde{O}(\sqrt{m})$ passes.

In Line 14, by Lemma 10.2, it can be done in one matrix multiplication with
$A_{3}$, which can be done in one pass and $O(n)$ space.

Outputting set $S$ requires $|S|$ space. ∎

######  Lemma 10.5 (Output slack stream).

Suppose there is an oracle $\mathcal{I}$ running in $f(n)$ space that can
output $b_{i}$ given any $i\in[m]$. Given a feasible $x\in\mathbb{R}^{n}$, we
can output a stream of all edges $e_{i}=(u,v)$ together with its slack
$s_{i}(x)=a_{i}^{\top}x-b_{i}$ in one pass within $O(n)+f(n)$ space.

###### Proof.

Without loss of generality, we assume every edge has its unique id $i\in[m]$
and we receive this id in the stream131313Actually we can remove this
assumption by observing that our oracle Algorithm 2 also works for
$i\in[n^{2}]$, in which we give every edge $e=(u,v)$ an index
$\mathrm{id}(e)=(u\cdot n+v)\in[n^{2}]$.. We read one pass of all edges. Upon
receiving edge $e_{i}=(u,v)$, we output $x_{u}-x_{v}-b_{i}$ where we use
$\mathcal{I}$ to calculate $b_{i}$. ∎

###  10.4 Building blocks

######  Lemma 10.6 (Hessian is SDDM matrix).

Let the graph $G=(V,E)$ be the input of Algorithm 5. Let $n=|V|$, $m=|E|$. Let
$A_{3}\in\mathbb{R}^{(m+n)\times n}$, $b_{3}\in\mathbb{R}^{m+n}$ be defined as
in Line 5 (Algorithm 5).

Let $x\in\mathbb{R}^{n}$ be any feasible point, i.e. $A_{3}x>b_{3}$. Let
$s(x)\in\mathbb{R}^{m+n}$ and $H(x)\in\mathbb{R}^{n\times n}$ be the slack and
hessian defined w.r.t. $A_{3}$ and $b_{3}$ (Definition 5.4). Then
$H(x)\in\mathbb{R}^{n\times n}$ can be written as

$\displaystyle H(x)=L(x)+D(x),$

where $L(x)\in\mathbb{R}^{n\times n}$ is the Laplacian matrix of graph $G$
with edge weight $\\{s_{1}(x)^{-2},\cdots,s_{m}(x)^{-2}\\}$,
$D(x)\in\mathbb{R}^{n\times n}$ is the diagonal matrix with each diagonal
entry $D(x)_{i,i}=s_{i+m}(x)^{-2}$.

###### Proof.

Let $A_{1}\in\mathbb{R}^{m\times n}$ be signed edge-vertex incident matrix of
$G$ with direction $V_{L}$ to $V_{R}$, then $A_{3}\cdot x\geq b_{3}$ can be
written as

$\displaystyle\begin{bmatrix}A_{1}\\\ I_{V_{L}}-I_{V_{R}}\end{bmatrix}\cdot
x\geq\begin{bmatrix}b_{1}\\\ 0\end{bmatrix},$

where $I_{V_{L}}\in\mathbb{R}^{n\times n}$ denotes the diagonal matrix with

$\displaystyle I_{i,i}=\begin{cases}1,&\mathrm{\leavevmode\nobreak\
if\leavevmode\nobreak\ }i\in V_{L};\\\ 0,&\mathrm{\leavevmode\nobreak\
otherwise}.\end{cases}$

and $I_{V_{R}}\in\mathbb{R}^{n\times n}$ denotes the diagonal matrix with

$\displaystyle I_{i,i}=\begin{cases}1,&\mathrm{\leavevmode\nobreak\
if\leavevmode\nobreak\ }i\in V_{R};\\\ 0,&\mathrm{\leavevmode\nobreak\
otherwise}.\end{cases}$

Denote $S\in\mathbb{R}^{m\times m}$ the diagonal matrix with each diagonal
entry $S_{i,i}:=s_{i}(x)$. Thus $H(x)$ can be written as

$\displaystyle H(x)=$ $\displaystyle\leavevmode\nobreak\
A_{1}^{\top}S^{-2}A_{1}+\sum_{i=m+1}^{m+n}e_{i}e_{i}^{\top}\cdot
s_{i}(x)^{-2}$ $\displaystyle=$ $\displaystyle\leavevmode\nobreak\ L(x)+D(x),$

where the first step is by definition of $H(x)$ (Definition 5.4), the second
step is by definition of $L(x)$ and $D(x)$. ∎

###  10.5 A minimum vertex cover solver

Note that Algorithm 5 is actually a high-accuracy fractional minimum vertex
cover solver for general graph (not necessarily bipartite graph), since we do
not use the bipartite property of matrix $A$ in IPM.

######  Theorem 10.7.

Let $G$ be a graph with $n$ vertices and $m$ edges. Consider the _fractional
minimum vertex cover_ problem Eq. (29) in which every edge $e$ needs to be
covered at least $b_{e}$ times. Let $x^{*}\in\mathbb{R}^{n}$ be the optimal
solution of Eq. (29). There exists a streaming algorithm (Algorithm 5) such
that for any $\delta>0$, it can output a feasible vertex cover
$x\in\mathbb{R}^{n}$ such that

$\mathbf{1}_{n}^{\top}x\leq\mathbf{1}_{n}^{\top}x^{*}+\delta$

in $\widetilde{O}(\sqrt{m})\cdot\log(1/\delta)$ passes and
$\widetilde{O}(n)\cdot\log(1/\delta)$ space with probability
$1-1/\operatorname{poly}(n)$.

###### Proof.

The proof follows directly from the proof of Lemma 10.1 and Lemma 10.4. ∎

As a byproduct, we obtain a fast semi-streaming algorithm for (exact,
integral) minimum vertex cover in bipartite graph.

######  Theorem 10.8.

Given a bipartite graph $G$ with $n$ vertices and $m$ edges, there exists a
streaming algorithm that computes a minimum vertex cover of $G$ in
$\widetilde{O}(\sqrt{m})$ passes and $\widetilde{O}(n)$ space with probability
$1-1/\operatorname{poly}(n)$.

###### Proof.

Because $G$ is bipartite, by Theorem 3.5, all extreme points of the polytope
of LP (29) are integral. Call Algorithm 5 to solve the perturbed LP in Lemma
10.2. Since LP (29) is feasible and bounded, the high-accuracy solution
obtained from IPM can be rounded to the optimal integral solution, which in
total takes $\widetilde{O}(\sqrt{m})$ passes and $\widetilde{O}(n)$ space with
probability $1-1/\operatorname{poly}(n)$. ∎

##  11 Combine

In this section, we give our main algorithm that combines all previous
subroutines to give our final Theorem 11.1. The proof consists of the
correctness part and the pass complexity part. We design the roadmap of this
section as follows. In Section 11.1, we present the main algorithm. In Section
11.2, we prove the pass complexity. The correctness is a bit complex. In
Section 11.4, we give preliminary knowledge on primal-to-dual transformation,
and then in section 11.5, combined with isolation lemma, we show how to
isolate a maximum weight matching from an optimal solution of the dual. In
Section 11.3, we conclude the correctness part.

######  Theorem 11.1 (Main theorem, formal version of Theorem 1.1).

Given a bipartite graph $G$ with $n$ vertices and $m$ edges, there exists a
streaming algorithm that computes a maximum weighted matching of $G$ in
$\widetilde{O}(\sqrt{m})$ passes and $\widetilde{O}(n)$ space with probability
$1-1/\operatorname{poly}(n)$.

###### Proof.

By combining Lemma 11.2 and Lemma 11.3. ∎

Statement | Section | Where we use | Comments
---|---|---|---
Theorem 11.1 | 11 | / | Final Theorem of this paper
Lemma 11.2 | 11.2 | Theorem 11.1 | Pass complexity
Lemma 11.3 | 11.3 | Theorem 11.1 | Correctness
Theorem 11.6 | 11.4 | Lemma 11.7 | Preliminary on proof of Lemma 11.3
Lemma 11.7 | 11.5 | Lemma 11.3 | Properties of Primal and Dual solutions
Table 5: The structure of Section 11.

###  11.1 Algorithms

Algorithm 6 is our main algorithm. The outer loop did $O(\log n)$ calls to
boost the success probability to $1-1/\operatorname{poly}(n)$. In each call,
it first prepares the isolation oracle from Section 4 and then passes it to
the minimum vertex cover solver from Section 10. The solver will return an
answer successfully with probability at least $1/2$. We choose the maximum
weighted matching over all $O(\log n)$ tries as the final answer.

Algorithm 6 $\textsc{Main}(G=(V_{L},V_{R},E,w))$

1:procedure Main($G=(V_{L},V_{R},E,w)$) $\triangleright$ Theorem 11.1

2: $\triangleright$ $w$ denotes edge weights that are integers.

3: $n\leftarrow|V_{L}|+|V_{R}|$, $m\leftarrow|E|$

4: $M\leftarrow\emptyset$ $\triangleright$ $M$ is maximum matching

5: for $i=1\to O(\log n)$ do

6: $\overline{b}\leftarrow\textsc{Isolation}(m,n^{n})$ $\triangleright$
$\overline{b}\in\operatorname*{\mathbb{Z}}^{m}$, Algorithm 1, Lemma 4.2

7: $b\leftarrow\overline{b}+n^{10}\cdot w$

8: $c\leftarrow\mathbf{1}_{V_{L}}-\mathbf{1}_{V_{R}}$

9: $S\leftarrow\textsc{MinimumVertexCover}(G,b,c,n,m)$ $\triangleright$
$S\subseteq[m]$, Algorithm 5

10: if $|S|\leq n$ then

11: Let $M^{\prime}$ be maximum weighted matching found in edge set $S$

12: if $w(M^{\prime})>w(M)$ then $\triangleright$ Update maximum weighted
matching

13: $M\leftarrow M^{\prime}$

14: end if

15: end if

16: end for

17: return $M$

18:end procedure

###  11.2 Pass complexity

The goal of this section is to bound the pass number (Lemma 11.2).

######  Lemma 11.2 (Pass complexity).

Given a bipartite graph with $n$ vertices and $m$ edges, Algorithm 6 can be
implemented such that it runs in $\widetilde{O}(\sqrt{m})$ passes in streaming
model in $\widetilde{O}(n)$ space.

###### Proof.

In Line 6, 7, we actually do not explicitly calculate $\overline{b}$ and $b$
and stored them in memory. We instead use an oracle $\mathcal{I}$ stated in
Lemma 4.1 that gives $b_{i}$ bits by bits. So this step does not cost space.

In Line 8, we calculate and store $c\in\operatorname*{\mathbb{Z}}^{n}$ in
$\widetilde{O}(n)$ space.

In Line 9, we call MinimuVertexCover. In order to use Lemma 10.4, we need to
prove the following properties.

1\. By Lemma 4.2, $\|\overline{b}\|_{\infty}\leq n^{7}$ so $\|b\|_{\infty}\leq
n^{10}$. And $\|c\|_{\infty}=1$. So both $\|b\|_{\infty}$ and $\|c\|_{\infty}$
is polynomially bounded.

2\. By Lemma 4.1, there is an oracle $\mathcal{I}$ that uses
$O(\log(n^{n})+\log(m))=\widetilde{O}(n)$ space that can output
$\overline{b}_{i}$. Since each edge $e_{i}$ comes with its weight $w_{i}$ in
the stream, we can output $b_{i}=\overline{b}_{i}+n^{10}w_{i}$ when given
$i\in[m]$.

3\. We can always assume $|S|\leq n$ since otherwise we will never enter Line
10.

By applying Lemma 10.4, this call can be done in $\widetilde{O}(n)$ space and