Given an undirected graph \(G = (V, E, w)\), a Gomory-Hu tree \(T\) (Gomory and Hu, 1961) is a tree on \(V\) that preserves all-pairs mincuts of \(G\) exactly.

We give a simple algorithm for the dynamic approximate All-Pairs Shortest Paths (APSP) problem. Given a graph G = (V, E, l) with polynomially bounded edge lengths, our data structure processes |E| edge insertions and deletions in total time |E|^{1+o(1)} and provides query access to |E|^o(1)-approximate distances in time Õ(1) per query. We produce a data structure that mimics Thorup-Zwick distance oracles [Thorup and Zwick, 2005], but is dynamic and deterministic. Our algorithm selects a small number of pivot vertices. Then, for every other vertex, it reduces distance computation to maintaining distances to a small neighborhood around that vertex and to the nearest pivot. We maintain distances between pivots efficiently by representing them in a smaller graph and recursing. We maintain these smaller graphs by (a) reducing vertex count using the dynamic distance-preserving core graphs of Kyng-Meierhans-Probst Gutenberg [Kyng et al., 2024] in a black-box manner and (b) reducing edge-count using a dynamic spanner akin to Chen-Kyng-Liu-Meierhans-Probst Gutenberg [Chen et al., 2024]. Our dynamic spanner internally uses an APSP data structure. Choosing a large enough size reduction factor in the first step allows us to simultaneously bootstrap a spanner and a dynamic APSP data structure. Notably, our approach does not need expander graphs, an otherwise ubiquitous tool in derandomization.

Given an undirected, weighted graph $G=(V, E, w)$, a Gomory-Hu tree or cut tree (Gomory and Hu, 1961) is a tree T over the vertex set V such that for every pair of vertices $s, t \in V$, the ($s, t$) min-cut in T is also an ($s, t$) min-cut in G and has the same value. In this article, we give the first deterministic almost-linear-time algorithm for constructing a Gomory-Hu tree. Our algorithm runs in $m^{1+o(1)}$-time, where m denotes the number of edges in the input graph G; this is clearly optimal up to the $m^{o(1)}$ term in the running time. Prior to our work, the best deterministic algorithm for this problem dated back to the original algorithm of Gomory and Hu that runs in $n m^{1+o(1)}$ time using current maxflow algorithms. In fact, our algorithm is also the first almost-linear-time deterministic algorithm for even simpler problems, such as finding the k-edge-connected components of a graph. Our new result hinges on two separate and novel components that each introduce a distinct set of de-randomization tools of independent interest: - a deterministic reduction from the all-pairs min-cuts problem to the single-source min-cuts problem incurring only sub-polynomial overhead, and - a deterministic almost-linear time algorithm for the singlesource min-cuts problem.

We present a new and surprisingly simple analysis of random-shift decompositions-originally proposed by Miller, Peng, and Xu [SPAA’13]: We show that decompositions for exponentially growing scales $D=2^{0}, 2^{1}, \ldots, 2^{\log _{2}(\operatorname{diam}(G))}$, have a tight constant trade-off between distance-to-center and separation probability on average across the distance scales - opposed to a necessary $\Omega(\log n)$ trade-off for a single scale. This almost immediately yields a way to compute a tree T for graph G that preserves all graph distances with expected $O(\log n)$-stretch. This gives an alternative proof that obtains tight approximation bounds of the seminal result by Fakcharoenphol, Rao, and Talwar [STOC’03] matching the $\Omega(\log n)$ lower bound by Bartal [FOCS’96]. Our insights can also be used to refine the analysis of a simple $\ell_{1}$-oblivious routing proposed in [FOCS’22], yielding a tight $O(\log n)$ competitive ratio. Our algorithms for constructing tree embeddings and $\ell_{1}$ oblivious routings can be implemented in the sequential, parallel, and distributed settings with optimal work, depth, and rounds, up to polylogarithmic factors. Previously, fast algorithms with tight guarantees were not known for tree embeddings in parallel and distributed settings, and for $\ell_{1}$-oblivious routings, not even a fast sequential algorithm was known.

We propose a randomized multiplicative weight update (MWU) algorithm for 𝓁_{∞} regression that runs in Õ(n^{2+1/22.5} poly(1/ε)) time when ω = 2+o(1), improving upon the previous best Õ(n^{2+1/18} polylog(1/ε)) runtime in the low-accuracy regime. Our algorithm combines state-of-the-art inverse maintenance data structures with acceleration. In order to do so, we propose a novel acceleration scheme for MWU that exhibits stability and robustness, which are required for the efficient implementations of the inverse maintenance data structures. We also design a faster deterministic MWU algorithm that runs in Õ(n^{2+1/12}poly(1/ε)) time when ω = 2+o(1), improving upon the previous best Õ(n^{2+1/6} poly log(1/ε)) runtime in the low-accuracy regime. We achieve this by showing a novel stability result that goes beyond previously known works based on interior point methods (IPMs). Our work is the first to use acceleration and inverse maintenance together efficiently, finally making the two most important building blocks of modern structured convex optimization compatible.

We give a simple algorithm for maintaining a $n^{o(1)}$-approximate spanner $H$ of a graph $G$ with $n$ vertices as $G$ receives edge updates by reduction to the dynamic All-Pairs Shortest Paths (APSP) problem. Given an initially empty graph $G$, our algorithm processes $m$ insertions and $n$ deletions in total time $m^{1 + o(1)}$ and maintains an initially empty spanner $H$ with total recourse $n^{1 + o(1)}$. When the number of insertions is much larger than the number of deletions, this notably yields recourse sub-linear in the total number of updates. Our algorithm only has a single $O(\log n)$ factor overhead in runtime and approximation compared to the underlying APSP data structure. Therefore, future improvements for APSP will directly yield an improved dynamic spanner.

A partially dynamic graph is a graph that undergoes edge insertions or deletions, but not both. In this talk, I present a unifying framework that yields the first almost-optimal, almost-linear time algorithms for many well-studied problems on partially dynamic graphs [Chen-Kyng-Liu-Meierhans-Probst-Gutenberg, STOC'24; Brand-Chen-Kyng-Liu-Meierhans-Probst Gutenberg-Sachdevea, FOCS'24]. These problems include cycle detection, strongly connected components, s-t distances, transshipment, bipartite matching, maximum flow, and minimum-cost flow. We achieve this unification by solving the partially dynamic threshold minimum-cost flow problem. We solve these problems by combining a partially dynamic L1 interior point method (Brand-Liu-Sidford STOC'23) with powerful new data structures that solve fully-dynamic APSP and min-cut with sub-polynomial approximation quality and sub-polynomial update and query time.

We present an algorithm that computes exact maximum flows and minimum-cost flows on directed graphs with m edges and polynomially bounded integral demands, costs, and capacities in \(m^{1+o(1)}\) time. Our algorithm builds the flow through a sequence of \(m^{1+o(1)}\) approximate undirected minimum-ratio cycles, each of which is computed and processed in amortized \(m^{o(1)}\) time using a new dynamic graph data structure. Our framework extends to algorithms running in \(m^{1+o(1)}\) time for computing flows that minimize general edge-separable convex functions to high accuracy. This gives almost-linear time algorithms for several problems including entropy-regularized optimal transport, matrix scaling, p -norm flows, and p -norm isotonic regression on arbitrary directed acyclic graphs.

We give the first almost-linear total time algorithm for deciding if a flow of cost at most$F$still exists in a directed graph, with edge costs and capacities, undergoing decremental updates, i.e., edge deletions, capacity decreases, and cost increases. This implies almost-linear time algorithms for approximating the minimum-cost flow value and s-t distance on such decremental graphs. Our framework additionally allows us to maintain decremental strongly connected components in almost-linear time deterministically. These algorithms also improve over the current best known runtimes for statically computing minimum-cost flow, in both the randomized and deterministic settings. We obtain our algorithms by taking the dual perspective, which yields cut-based algorithms. More precisely, our algorithm computes the flow via a sequence of$m^{1+o(1)}$-dynamic min-ratio cut problems, the dual analog of the dynamic min-ratio cycle problem that underlies recent fast algorithms for minimum-cost flow. Our main technical contribution is a new data structure that returns an approximately optimal min-ratio cut in amortized$m^{o(1)}$time by maintaining a tree-cut sparsifier. This is achieved by devising a new algorithm to maintain the dynamic expander hierarchy of [$\text{Goranci-Racke-}$SaranurakTan, SODA 2021] that also works in capacitated graphs. All our algorithms are deterministc, though they can be sped up further using randomized techniques while still working against an adaptive adversary.

In this paper, we investigate the question of whether the electrical flow routing is a good oblivious routing scheme on an m-edge graph G = (V, E) that is a Φ-expander, i.e. where |∂ S| ≥ Φ ⋅ vol(S) for every S ⊆ V, vol(S) ≤ vol(V)/2. Beyond its simplicity and structural importance, this question is well-motivated by the current state-of-the-art of fast algorithms for 𝓁_∞ oblivious routings that reduce to the expander-case which is in turn solved by electrical flow routing. Our main result proves that the electrical routing is an O(Φ^{-1} log m)-competitive oblivious routing in the 𝓁₁- and 𝓁_∞-norms. We further observe that the oblivious routing is O(log² m)-competitive in the 𝓁₂-norm and, in fact, O(log m)-competitive if 𝓁₂-localization is O(log m) which is widely believed. Using these three upper bounds, we can smoothly interpolate to obtain upper bounds for every p ∈ [2, ∞] and q given by 1/p + 1/q = 1. Assuming 𝓁₂-localization in O(log m), we obtain that in 𝓁_p and 𝓁_q, the electrical oblivious routing is O(Φ^{-(1-2/p)}log m) competitive. Using the currently known result for 𝓁₂-localization, this ratio deteriorates by at most a sublogarithmic factor for every p, q ≠ 2. We complement our upper bounds with lower bounds that show that the electrical routing for any such p and q is Ω(Φ^{-(1-2/p)} log m)-competitive. This renders our results in 𝓁₁ and 𝓁_∞ unconditionally tight up to constants, and the result in any 𝓁_p- and 𝓁_q-norm to be tight in case of 𝓁₂-localization in O(log m).

We provide an algorithm which, with high probability, maintains a (1 — ɛ)-approximate maximum flow on an undirected graph undergoing m-edge additions in amortized mo(1)ɛ-3 time per update. To obtain this result, we provide a more general algorithm that solves what we call the incremental, thresholded, p-norm flow problem that asks to determine the first edge-insertion in an undirected graph that causes the minimum ℓp-norm flow to decrease below a given threshold in value. Since we solve this thresholded problem, our data structure succeeds against an adaptive adversary that can only see the data structure's output. Furthermore, since our algorithm holds for p = 2, we obtain improved algorithms for dynamically maintaining the effective resistance between a pair of vertices in an undirected graph undergoing edge insertions.

In a breakthrough result, Spielman and Teng (2004) developed a nearly-linear time solver for Laplacian linear equations, i.e. equations where the coefficient matrix is symmetric with non-negative diagonals and zero row sums. Since the development of the Spielman-Teng solver, there has been substantial progress, simplifying and improving their result, but obtaining a fast practical, parallel Laplacian solver remains an open problem.

We give the first almost-linear time algorithms for several problems in incremental graphs including cycle detection, strongly connected component maintenance, s-t shortest path, maximum flow, and minimum-cost flow. To solve these problems, we give a deterministic data structure that returns a mo(1)-approximate minimum-ratio cycle in fully dynamic graphs in amortized mo(1) time per update. Combining this with the interior point method framework of Brand-Liu-Sidford (STOC 2023) gives the first almost-linear time algorithm for deciding the first update in an incremental graph after which the cost of the minimum-cost flow attains value at most some given threshold F. By rather direct reductions to minimum-cost flow, we are then able to solve the problems in incremental graphs mentioned above.

We present a general toolbox, based on new vertex sparsifiers, for designing data structures to maintain shortest paths in graphs undergoing edge insertions and/or deletions. In particular, we obtain the following results:

The \(\ell _p\) -norm regression problem is a classic problem in optimization with wide ranging applications in machine learning and theoretical computer science. The goal is to compute \(\boldsymbol {\mathit {x}}^{\star } =\arg \min _{\boldsymbol {\mathit {A}}\boldsymbol {\mathit {x}}=\boldsymbol {\mathit {b}}}\Vert \boldsymbol {\mathit {x}}\Vert _p^p\) , where \(\boldsymbol {\mathit {x}}^{\star }\in \mathbb {R}^n,\boldsymbol {\mathit {A}}\in \mathbb {R}^{d\times n},\boldsymbol {\mathit {b}}\in \mathbb {R}^d\) and \(d\le n\) . Efficient high-accuracy algorithms for the problem have been challenging both in theory and practice and the state-of-the-art algorithms require \(poly(p)\cdot n^{\frac{1}{2}-\frac{1}{p}}\) linear system solves for \(p\ge 2\) . In this article, we provide new algorithms for \(\ell _p\) -regression (and a more general formulation of the problem) that obtain a high-accuracy solution in \(O(p n^{ {(p-2)}{(3p-2)}})\) linear system solves. We further propose a new inverse maintenance procedure that speeds-up our algorithm to \(\widetilde{O}(n^{\omega })\) total runtime, where \(O(n^{\omega })\) denotes the running time for multiplying \(n \times n\) matrices. Additionally, we give the first Iteratively Reweighted Least Squares (IRLS) algorithm that is guaranteed to converge to an optimum in a few iterations. Our IRLS algorithm has shown exceptional practical performance, beating the currently available implementations in MATLAB/CVX by 10–50×.

We give a deterministic $m^{1+o(1)}$ time algorithm that computes exact maximum flows and minimum-cost flows on directed graphs with m edges and polynomially bounded integral demands, costs, and capacities. As a consequence, we obtain the first running time improvement for deterministic algorithms that compute maximum-flow in graphs with polynomial bounded capacities since the work of Goldberg-Rao [J.ACM ’98].Our algorithm builds on the framework of Chen-Kyng-Liu-Peng-Gutenberg-Sachdeva [FOCS ’22] that computes an optimal flow by computing a sequence of $m^{1+o(1)}$-approximate undirected minimum-ratio cycles. We develop a deterministic dynamic graph data-structure to compute such a sequence of minimum-ratio cycles in an amortized $m^{o(1)}$ time per edge update. Our key technical contributions are deterministic analogues of the vertex sparsification and edge sparsification components of the data-structure from Chen et al. For the vertex sparsification component, we give a method to avoid the randomness in Chen et al. which involved sampling random trees to recurse on. For the edge sparsification component, we design a deterministic algorithm that maintains an embedding of a dynamic graph into a sparse spanner. We also show how our dynamic spanner can be applied to give a deterministic data structure that maintains a fully dynamic low-stretch spanning tree on graphs with polynomially bounded edge lengths, with subpolynomial average stretch and subpolynomial amortized time per edge update.

In this article, we show that the algorithm of maintaining expander decompositions in graphs undergoing edge deletions directly by removing sparse cuts repeatedly can be made efficient. Formally, for an m-edge undirected graph G, we say a cut is ϕ-sparse if . A ϕ-expander decomposition of G is a partition of V into sets X1,X2,…, Xk such that each cluster G[X1] contains no ϕ-sparse cut (meaning it is a ϕ-expander) with Õ(ϕm) edges crossing between clusters. A natural way to compute a ϕ-expander decomposition is to decompose clusters by ϕ-sparse cuts until no such cut is contained in any cluster. We show that even in graphs undergoing edge deletions, a slight relaxation of this meta-algorithm can be implemented efficiently with amortized update time mo(1)/ϕ2. Our approach naturally extends to maintaining directed ϕ-expander decompositions and ϕ-expander hierarchies and thus gives a unifying framework while having simpler proofs than previous state-of-the-art work. In all settings, our algorithm matches the run-times of previous algorithms up to subpolynomial factors. Moreover, our algorithm provides stronger guarantees for ϕ-expander decompositions. For example, for graphs undergoing edge deletions, our approach is the first to maintain a dynamic expander decomposition where each updated decomposition is a refinement of the previous decomposition, and our approach is the first to guarantee a sublinear ϕm1+ο(1) bound on the total number of edges that cross between clusters across the entire sequence of dynamic updates. Our techniques also give by far the simplest, deterministic algorithms for maintaining Strongly-Connected Components (SCCs) in directed graphs undergoing edge deletions, and for maintaining connectivity in undirected fully-dynamic graphs, both matching the current state-of-the art run-times up to subpolynomial factors. * The full version of the paper can be accessed at https://arxiv.org/abs/2204.02519

We prove lower bounds for higher-order methods in smooth non-convex finite-sum optimization. Our contribution is threefold: We first show that a deterministic algorithm cannot profit from the finite-sum structure of the objective and that simulating a pth-order regularized method on the whole function by constructing exact gradient information is optimal up to constant factors. We further show lower bounds for randomized algorithms and compare them with the best-known upper bounds. To address some gaps between the bounds, we propose a new second-order smoothness assumption that can be seen as an analogue of the first-order mean-squared smoothness assumption. We prove that it is sufficient to ensure state-of-the-art convergence guarantees while allowing for a sharper lower bound.

We give an algorithm that computes exact maximum flows and minimum-cost flows on directed graphs with m edges and polynomially bounded integral demands, costs, and capacities in $m^{1+o(1)}$ time. Our algorithm builds the flow through a sequence of $m^{1+o(1)}$ approximate undirected minimum-ratio cycles, each of which is computed and processed in amortized $m^{o(1)}$ time using a new dynamic graph data structure. Our framework extends to algorithms running in $m^{1+o(1)}$ time for computing flows that minimize general edge-separable convex functions to high accuracy. This gives almost-linear time algorithms for several problems including entropy-regularized optimal transport, matrix scaling, p-norm flows, and p-norm isotonic regression on arbitrary directed acyclic graphs.

In this article, we present the first deterministic directed Laplacian L systems solver that runs in time almost-linear in the number of non-zero entries of L. Previous reductions imply the first deterministic almost-linear time algorithms for computing various fundamental quantities on directed graphs including stationary distributions, personalized PageRank, hitting times and escape probabilities.We obtain these results by introducing partial symmetrization, a new technique that makes the Laplacian of an Eulerian directed graph “less directed” in a useful sense, which may be of independent interest. The usefulness of this technique comes from two key observations: Firstly, the partially symmetrized Laplacian preconditions the original Eulerian Laplacian well in Richardson iteration, enabling us to construct a solver for the original matrix from a solver for the partially symmetrized one. Secondly, the undirected structure in the partially symmetrized Laplacian makes it possible to sparsify the matrix very crudely, i.e. with large spectral error, and still show that Richardson iterations convergence when using the sparsified matrix as a preconditioner. This allows us to develop deterministic sparsification tools for the partially symmetrized Laplacian.Together with previous reductions from directed Laplacians to Eulerian Laplacians, our technique results in the first deterministic almost-linear time algorithm for solving linear equations in directed Laplacians. To emphasize the generality of our new technique, we show that two prominent existing (randomized) frameworks for solving linear equations in Eulerian Laplacians can be derandomized in this way: the squaring-based framework of [1] and the sparsified Cholesky-based framework of [2].

We give a nearly-linear time reduction that encodes any linear program as a 2-commodity flow problem with only a small blow-up in size. Under mild assumptions similar to those employed by modern fast solvers for linear programs, our reduction causes only a polylogarithmic multiplicative increase in the size of the program and runs in nearly-linear time. Our reduction applies to high-accuracy approximation algorithms and exact algorithms. Given an approximate solution to the 2-commodity flow problem, we can extract a solution to the linear program in linear time with only a polynomial factor increase in the error. This implies that any algorithm that solves the 2-commodity flow problem can solve linear programs in essentially the same time. Given a directed graph with edge capacities and two source-sink pairs, the goal of the 2-commodity flow problem is to maximize the sum of the flows routed between the two source-sink pairs subject to edge capacities and flow conservation. A 2-commodity flow can be directly written as a linear program, and thus we establish a nearly-tight equivalence between these two classes of problems. Our proof follows the outline of Itai's polynomial-time reduction of a linear program to a 2-commodity flow problem (JACM'78). Itai's reduction shows that exactly solving 2-commodity flow and exactly solving linear programming are polynomial-time equivalent. We improve Itai's reduction to nearly preserve the problem representation size in each step. In addition, we establish an error bound for approximately solving each intermediate problem in the reduction, and show that the accumulated error is polynomially bounded. We remark that our reduction does not run in strongly polynomial time and that it is open whether 2-commodity flow and linear programming are equivalent in strongly polynomial time.

We study linear equations in combinatorial Laplacians of k-dimensional simplicial complexes (k-complexes), a natural generalization of graph Laplacians. Combinatorial Laplacians play a crucial role in homology and are a central tool in topology. Beyond this, they have various applications in data analysis and physical modeling problems. It is known that nearly-linear time solvers exist for graph Laplacians. However, nearly-linear time solvers for combinatorial Laplacians are only known for restricted classes of complexes. This paper shows that linear equations in combinatorial Laplacians of 2-complexes are as hard to solve as general linear equations. More precisely, for any constant c ≥ 1, if we can solve linear equations in combinatorial Laplacians of 2-complexes up to high accuracy in time Õ((# of nonzero coefficients)^c), then we can solve general linear equations with polynomially bounded integer coefficients and condition numbers up to high accuracy in time Õ((# of nonzero coefficients)^c). We prove this by a nearly-linear time reduction from general linear equations to combinatorial Laplacians of 2-complexes. Our reduction preserves the sparsity of the problem instances up to poly-logarithmic factors.

We present new scalar and matrix Chernoff-style concentration bounds for a broad class of probability distributions over the binary hypercube {0, 1}n. Motivated by recent tools developed for the study of mixing times of Markov chains on discrete distributions, we say that a distribution is ℓ∞-independent when the infinity norm of its influence matrix is bounded by a constant. We show that any distribution which is ℓ∞-infinity independent satisfies a matrix Chernoff bound that matches the matrix Chernoff bound for independent random variables due to Tropp. Our matrix Chernoff bound is a broad generalization and strengthening of the matrix Chernoff bound of Kyng and Song (FOCS'18). Using our bound, we can conclude as a corollary that a union of O(log |V|) random spanning trees gives a spectral graph sparsifier of a graph with |V| vertices with high probability matching results for independent edge sampling, and matching lower bounds from Kyng and Song.

Given a directed, weighted graph G = (V, E) undergoing edge insertions, the incremental single-source shortest paths (SSSP) problem asks for the maintenance of approximate distances from a dedicated source s while optimizing the total time required to process the insertion sequence of m edges. Recently, Gutenberg, Williams and Wein [STOC'20] introduced a deterministic Õ(n2) algorithm for this problem, achieving near linear time for very dense graphs. For sparse graphs, Chechik and Zhang [SODA'21] recently presented a deterministic Õ(m5/3) algorithm, and an adaptive randomized algorithm with run-time . This algorithm is remarkable for two reasons: 1) in very spare graphs it reaches the directed hopset barrier of Ω(n3/2) that applied to all previous approaches for partially-dynamic SSSP [STOC'14, SODA'20, FOCS'20] and 2) it does not resort to a directed hopset technique itself. In this article we introduce propagation synchronization, a new technique for controlling the error build-up on paths throughout batches of insertions. This leads us to a significant improvement of the approach in [SODA'21] yielding a deterministic O(m3/2) algorithm for the problem. By a very careful combination of our new technique with the sampling approach from [SODA'21], we further obtain an adaptive randomized algorithm with total update time Õ(m4/3). This is the first partially-dynamic SSSP algorithm in sparse graphs to bypass the notorious directed hopset barrier which is often seen as the fundamental challenge towards achieving truly near-linear time algorithms.

We give almost-linear-time algorithms for constructing sparsifiers with $n\ poly(\log n)$ edges that approximately preserve weighted $(\ell^{2}_2 + \ell^{p}_p)$ flow or voltage objectives on graphs. For flow objectives, this is the first sparsifier construction for such mixed objectives beyond unit $\ell_p$ weights, and is based on expander decompositions. For voltage objectives, we give the first sparsifier construction for these objectives, which we build using graph spanners and leverage score sampling. Together with the iterative refinement framework of [Adil et al, SODA 2019], and a new multiplicative-weights based constant-approximation algorithm for mixed-objective flows or voltages, we show how to find $(1+2^{-\text{poly}(\log n)})$ approximations for weighted $\ell_p$-norm minimizing flows or voltages in $p(m^{1+o(1)} + n^{4/3 + o(1)})$ time for $p=ω(1),$ which is almost-linear for graphs that are slightly dense ($m \ge n^{4/3 + o(1)}$).

We study the complexity of approximately solving packing linear programs. In the Real RAM model, it is known how to solve packing LPs with N non-zeros in time Õ(N/ϵ). We investigate whether the ϵ dependence in the running time can be improved. Our first main result relates the difficulty of this problem to hardness assumptions for solving dense linear equations. We show that, in the Real RAM model, unless linear equations in matrices n × n with condition number O(n10) can be solved to ϵ accuracy faster than Õ(n2.01 log(1/ϵ)), no algorithm (1−ϵ)-approximately solves a O(n)×O(n) packing LPs (where N = O(n2)) in time Õ(n2ϵ−0.0003). It would be surprising to solve linear equations in the Real RAM model this fast, as we currently cannot solve them faster than Õ(nω), where ω denotes the exponent in the running time for matrix multiplication in the Real RAM model (and equivalently matrix inversion). The current best bound on this exponent is roughly ω ≤ 2.372. Note, however, that a fast solver for linear equations does not directly imply faster matrix multiplication. But, our reduction shows that if fast and accurate packing LP solvers exist, then either linear equations can be solved much faster than matrix multiplication or the matrix multiplication constant is very close to 2. Instantiating the same reduction with different parameters, we show that unless linear equations in matrices with condition number O(n1.5) can be solved to ϵ accuracy faster than Õ(n2.372 log(1/ϵ)), no algorithm (1 – ϵ)-approximately solves packing LPs in time Õ(n2ϵ−0.067). Thus smaller improvements in the exponent for ϵ in the running time of Packing LP solvers also imply improvements in the current state-of-the-art for solving linear equations. Our second main result relates the difficulty of approximately solving packing linear programs to hardness assumptions for solving sparse linear equations: In the Real RAM model, unless well-conditioned sparse systems of linear equations can be solved faster than Õ((no. non-zeros of matrix) ), no algorithm (1 – ϵ)-approximately solves packing LPs with N non-zeros in time Õ(Nϵ−0.165). This running time of Õ((no. non-zeros of matrix) ) is obtained by the classical Conjugate Gradient algorithm by a standard analysis. Our reduction implies that if sufficiently good packing LP solvers exist, then this long-standing best-known bound on the running time for solving well-conditioned systems of linear equations is sub-optimal1. While we prove results in the Real RAM model, our condition number assumptions ensure that our results can be translated to fixed point arithmetic with (log n)O(1) bits per number.

We show that if the nearly linear time solvers for Laplacian matrices and their generalizations can be extended to solve just slightly larger families of linear systems, then they can be used to quickly solve all systems of linear equations over the reals. This result can be viewed either positively or negatively: either nearly linear time algorithms can be developed for solving all systems of linear equations over the reals, or progress on the families that can be solved in nearly linear time will soon halt.

We give improved algorithms for the ℓp-regression problem, minx ‖x‖p such that Ax = b, for all p ∊ (1, 2) ∪ (2, ∞). Our algorithms obtain a high accuracy solution in iterations, where each iteration requires solving an m × m linear system, with m being the dimension of the ambient space. Incorporating a procedure for maintaining an approximate inverse of the linear systems that we need to solve at each iteration, we give algorithms for solving ℓp-regression to 1/poly(n) accuracy that runs in time Õp(mmax{ω, 7/3}), where ω is the matrix multiplication constant. For the current best value of ω > 2.37, this means that we can solve ℓp regression as fast as ℓ2 regression, for all constant p bounded away from 1. Our algorithms can be combined with nearly-linear time solvers for linear systems in graph Laplacians to give minimum ℓp-norm flow / voltage solutions to 1/poly(n) accuracy on an undirected graph with m edges in time. For sparse graphs and for matrices with similar dimensions, our iteration counts and running times improve upon the p-norm regression algorithm by [Bubeck-Cohen-Lee-Li STOC'18], as well as general purpose convex optimization algorithms. At the core of our algorithms is an iterative refinement scheme for ℓp-norms, using the quadratically-smoothed ℓp-norms introduced in the work of Bubeck et al. Formally, given an initial solution, we construct a problem that seeks to minimize a quadratically-smoothed ℓp norm over a subspace, such that a crude solution to this problem allows us to improve the initial solution by a constant factor, leading to algorithms with fast convergence.

We present algorithms for solving a large class of flow and regression problems on unit weighted graphs to (1 + 1 / poly(n)) accuracy in almost-linear time. These problems include ℓp-norm minimizing flow for p large (p ∈ [ω(1), o(log2/3 n) ]), and their duals, ℓp-norm semi-supervised learning for p close to 1.

In this paper, we show how to solve directed Laplacian systems in nearly-linear time. Given a linear system in an n × n Eulerian directed Laplacian with m nonzero entries, we show how to compute an ε-approximate solution in time O(m logO(1)(n) log (1/ε)). Through reductions from [Cohen et al. FOCS'16], this gives the first nearly-linear time algorithms for computing ε-approximate solutions to row or column diagonally dominant linear systems (including arbitrary directed Laplacians) and computing ε-approximations to various properties of random walks on directed graphs, including stationary distributions, personalized PageRank vectors, hitting times, and escape probabilities. These bounds improve upon the recent almost-linear algorithms of [Cohen et al. STOC'17], which gave an algorithm to solve Eulerian Laplacian systems in time O((m+n2O(√ log n log log n))logO(1)(n ε-1)). To achieve our results, we provide a structural result that we believe is of independent interest. We show that Eulerian Laplacians (and therefore the Laplacians of all strongly connected directed graphs) have sparse approximate LU-factorizations. That is, for every such directed Laplacian there are lower upper triangular matrices each with at most Õ(n) nonzero entries such that there product spectrally approximates the directed Laplacian in an appropriate norm. This claim can be viewed as an analog of recent work on sparse Cholesky factorizations of Laplacians of undirected graphs. We show how to construct such factorizations in nearly-linear time and prove that once constructed they yield nearly-linear time algorithms for solving directed Laplacian systems.

Strongly Rayleigh distributions are a class of negatively dependent distributions of binary-valued random variables [Borcea, Brändén, Liggett JAMS 09]. Recently, these distributions have played a crucial role in the analysis of algorithms for fundamental graph problems, e.g. Traveling Salesman Problem [Gharan, Saberi, Singh FOCS 11]. We prove a new matrix Chernoff bound for Strongly Rayleigh distributions. As an immediate application, we show that adding together the Laplacians of ε-2log2n random spanning trees gives an (1± ε) spectral sparsifiers of graph Laplacians with high probability. Thus, we positively answer an open question posted in [Baston, Spielman, Srivastava, Teng JACM 13]. Our number of spanning trees for spectral sparsifier matches the number of spanning trees required to obtain a cut sparsifier in [Fung, Hariharan, Harvey, Panigraphi STOC 11]. The previous best result was by naively applying a classical matrix Chernoff bound which requires ε-2n log n spanning trees. For the tree averaging procedure to agree with the original graph Laplacian in expectation, each edge of the tree should be reweighted by the inverse of the edge leverage score in the original graph. We also show that when using this reweighting of the edges, the Laplacian of single random tree is bounded above in the PSD order by the original graph Laplacian times a factor log n with high probability, i.e. LT≤ O(log n) LG. We show a lower bound that almost matches our last result, namely that in some graphs, with high probability, the random spanning tree is not bounded above in the spectral order by [log n/loglog n] times the original graph Laplacian. We also show a lower bound that in ε-2log n spanning trees are necessary to get a (1± ε) spectral sparsifier.

We present an asymptotically faster algorithm for solving linear systems in well-structured 3-dimensional truss stiffness matrices. These linear systems arise from linear elasticity problems, and can be viewed as extensions of graph Laplacians into higher dimensions. Faster solvers for the 2-D variants of such systems have been studied using generalizations of tools for solving graph Laplacians [Daitch-Spielman CSC’07, Shklarski-Toledo SIMAX’08].

We show that if the nearly-linear time solvers for Laplacian matrices and their generalizations can be extended to solve just slightly larger families of linear systems, then they can be used to quickly solve all systems of linear equations over the reals. This result can be viewed either positively or negatively: either we will develop nearly-linear time algorithms for solving all systems of linear equations over the reals, or progress on the families we can solve in nearly-linear time will soon halt.

We develop a framework for graph sparsification and sketching, based on a new tool, short cycle decomposition, which is a decomposition of an unweighted graph into an edge-disjoint collection of short cycles, plus a small number of extra edges. A simple observation shows that every graph $G$ on $n$ vertices with $m$ edges can be decomposed in $O(mn)$ time into cycles of length at most $2 \log n,$ and at most $2n$ extra edges. We give an $(m^{1+o(1)})$-time algorithm for constructing a short cycle decomposition, with cycles of length $n^{o(1)},$ and $n^{1+o(1)}$ extra edges. Both the existential and algorithmic variants of this decomposition enable us to make the following progress on several open problems in randomized graph algorithms: (1) We present an algorithm that runs in time $m^{1+o(1)}\varepsilon^{-1.5}$ and returns $(1\pm\varepsilon)$-approximations to effective resistances of all edges, improving over the previous best runtime of $\widetilde{{O}}(\min\{m\varepsilon^{-2}, n^{2} \varepsilon^{-1}\})$. This routine in turn gives an algorithm for approximating the determinant of a graph Laplacian up to a factor of $(1\pm \varepsilon)$ in $m^{1 + o(1)} + n^{\nicefrac{15}{8}+o(1)}\varepsilon^{-\nicefrac{7}{4}}$ time. (2) We show the existence of graphical spectral sketches with about $n\varepsilon^{-1}$ edges, and also give efficient algorithms to construct them. A graphical spectral sketch is a distribution over sparse graphs $H$ such that for a fixed vector ${\mathit{x}}$, we have ${{x}}^{\top} {L}_H {{x}} = (1\pm\varepsilon) {{x}}^{\top} {L}_G {{x}}$ and ${{x}}^{\top} {L}^{+}_H {{x}} = (1\pm\varepsilon) {{x}}^{\top} {L}^{+}_G {{x}}$ with high probability, where ${L}$ is the graph Laplacian and ${L}^{+}$ is its pseudoinverse. This implies the existence of resistance sparsifiers with about $n \varepsilon^{-1}$ edges that preserve the effective resistance between every pair of vertices up to $(1\pm\varepsilon)$. (3) By combining short cycle decompositions with known tools in graph sparsification, we show the existence of nearly linear sized degree-preserving spectral sparsifiers, as well as significantly sparser approximations of Eulerian directed graphs. The latter is critical to recent breakthroughs on faster algorithms for solving linear systems in directed Laplacians. The running time and output qualities of our spectral sketch and degree-preserving (directed) sparsification algorithms are limited by the efficiency of our routines for constructing short cycle decompositions. Improved algorithms for short cycle decompositions will lead to improvement in each of these algorithms.

We present an algorithm that, with high probability, generates a random spanning tree from an edge-weighted undirected graph in (n5/3 m1/3) time. The tree is sampled from a distribution where the probability of each tree is proportional to the product of its edge weights. This improves upon the previous best algorithm due to Colbourn et al. that runs in matrix multiplication time, O(nω). For the special case of unweighted graphs, this improves upon the best previously known running time of Õ(min{nω,m√n,m4/3}) for m ⪢ n7/4 (Colbourn et al. '96, Kelner-Madry '09, Madry et al. '15).

We show how to perform sparse approximate Gaussian elimination for Laplacian matrices. We present a simple, nearly linear time algorithm that approximates a Laplacian by the product of a sparse lower triangular matrix with its transpose. This gives the first nearly linear time solver for Laplacian systems that is based purely on random sampling, and does not use any graph theoretic constructions such as low-stretch trees, sparsifiers, or expanders. Our algorithm performs a subsampled Cholesky factorization, which we analyze using matrix martingales. As part of the analysis, we give a proof of a concentration inequality for matrix martingales where the differences are sums of conditionally independent variables.

We introduce the sparsified Cholesky and sparsified multigrid algorithms for solving systems of linear equations. These algorithms accelerate Gaussian elimination by sparsifying the nonzero matrix entries created by the elimination process. We use these new algorithms to derive the first nearly linear time algorithms for solving systems of equations in connection Laplacians---a generalization of Laplacian matrices that arise in many problems in image and signal processing. We also prove that every connection Laplacian has a linear sized approximate inverse. This is an LU factorization with a linear number of nonzero entries that is a strong approximation of the original matrix. Using such a factorization one can solve systems of equations in a connection Laplacian in linear time. Such a factorization was unknown even for ordinary graph Laplacians.

We develop fast algorithms for solving regression problems on graphs where one is given the value of a function at some vertices, and must find its smoothest possible extension to all vertices. The extension we compute is the absolutely minimal Lipschitz extension, and is the limit for large p of p-Laplacian regularization. We present an algorithm that computes a minimal Lipschitz extension in expected linear time, and an algorithm that computes an absolutely minimal Lipschitz extension in expected time \widetildeO (m n). The latter algorithm has variants that seem to run much faster in practice. These extensions are particularly amenable to regularization: we can perform l_0-regularization on the given values in polynomial time and l_1-regularization on the initial function values and on graph edge weights in time \widetildeO (m^3/2). Our definitions and algorithms naturally extend to directed graphs.

Given a directed acyclic graph $G,$ and a set of values $y$ on the vertices, the Isotonic Regression of $y$ is a vector $x$ that respects the partial order described by $G,$ and minimizes $\|x-y\|,$ for a specified norm. This paper gives improved algorithms for computing the Isotonic Regression for all weighted $\ell_{p}$-norms with rigorous performance guarantees. Our algorithms are quite practical, and their variants can be implemented to run fast in practice.

We show an algorithm for solving symmetric diagonally dominant (SDD) linear systems with m non-zero entries to a relative error of ε in O(m log1/2 n logc n log(1/ε)) time. Our approach follows the recursive preconditioning framework, which aims to reduce graphs to trees using iterative methods. We improve two key components of this framework: random sampling and tree embeddings. Both of these components are used in a variety of other algorithms, and our approach also extends to the dual problem of computing electrical flows.