Truly Subquadratic 3SUM and Truly Subcubic APSP via Triangles in Sparse Lopsided Graphs

Claude, an AI model developed by Anthropic, discovered the algorithm that refutes the 3SUM, APSP, and Exact Triangle hypotheses. The authors then worked to understand, simplify, strengthen, and extend the algorithm, derive additional consequences, and make the presentation accessible. See “Acknowledgments and Methodology” for how the result was found and shared with the authors. The authors take full responsibility for this paper.

Claude also verified this paper’s main results using the Lean 4 proof assistant with the Mathlib library.

Introduction

Since its inception, algorithms research has sought the fastest possible algorithm for every problem it considers. Powerful techniques have been developed, from variants of dynamic programming to surprising algebraic tools such as the Fast Fourier Transform [65] and fast matrix multiplication [147]. Yet for many problems, these techniques have not sufficed. Some examples include:

• 3SUM: Given $n$ numbers, decide whether three of them sum to 0. This is a classical problem with a long history, and it is central to computational geometry (see [92]). Despite many decades of research, the $O(n^2)$-time algorithm taught in algorithms classes has only been sped up by polylogarithmic factors [35, 94, 55].

• All-Pairs Shortest Paths (APSP): Given an edge-weighted graph on $n$ vertices with no negative cycles, compute the shortest-path distance between every pair of vertices. Classical algorithms such as Johnson’s [107] and Floyd and Warshall’s [85, 162] solve APSP in $O(n^3)$ time. After a long line of work [87, 150, 76, 97, 151, 152, 173, 52, 98, 53, 102], the fastest known algorithm [166] only shaves a subpolynomial, $2^{\Omega(\sqrt{\log n})}$ factor off the cubic running time.

• Directed Unweighted APSP: The version of APSP in $n$-vertex directed unweighted graphs. The fastest known algorithm, by Zwick [172], runs in $O(n^{2.5275})$ time with the current bounds for rectangular matrix multiplication [10], and beyond improvements to the matrix multiplication bounds, it has not been improved for more than two decades. If the exponent $\omega$ of matrix multiplication is 2, the best known running time for this problem [14, 172] is still only $\widetilde{O}(n^{2.5})$.

• Tree Edit Distance: Given two rooted, ordered, labeled trees with $n$ vertices, and costs for deleting, inserting, and relabeling vertices, find the cheapest sequence of such operations that transforms one tree into the other. Tree Edit Distance is used, for instance, to compare RNA secondary structures and XML documents. The fastest known algorithms for this problem [75, 133] run in cubic time, up to subpolynomial improvements. Only when all costs are equal to 1 are truly subcubic algorithms known [127, 74, 133].

• Online Set Disjointness: Preprocess sets $S_1,\ldots,S_N$ and $T_1,\ldots,T_N$ over a universe of size $D$, so that given $I$ and $J$ one can quickly decide whether $S_I \cap T_J = \emptyset$, or, in the counting version, return $\lvert S_I \cap T_J\rvert$. Essentially only two straightforward solutions are known: (1) to compute all $N^2$ answers in advance, which fast rectangular matrix multiplication does in $N^{2+o(1)}$ time when $D \le N^{0.321}$ [64, 161], and (2) to store the sets and compute each intersection at query time in $O(D)$ time. Set Disjointness data structures have been studied extensively [50, 116, 89, 155, 119], mostly for sets of small total size over a large universe.

• CNF-SAT: Given a Boolean formula $F$ in conjunctive normal form over $n$ variables and $O(n)$ clauses, is there an assignment to the variables that satisfies all clauses? The brute-force algorithm solves CNF-SAT in $2^n\operatorname{poly}(n)$ time. There are slightly faster algorithms when the width of the clauses is bounded by a constant independent of $n$ [144, 136, 99]. However, no $O(1.99999^n)$-time algorithm is known.

• Orthogonal Vectors (OV): Given $n$ Boolean vectors in $d = \omega(\log n)$ dimensions, are two of them orthogonal? The brute-force algorithm solves OV in $O(n^2d)$ time, and only subpolynomial improvements are known: the fastest algorithms [29, 66] run in $n^{2-1/O(\log(d/\log n))}$ time.

The apparent lack of progress on problems such as these led to new theories of hardness: NP-completeness [61, 110] for problems with no known polynomial-time algorithm, and fine-grained complexity (FGC) [154, 45] for problems whose best known algorithms were stuck at essentially the straightforward running time. Both theories relate problems via reductions. FGC focuses on improvements in the exponent: a fine-grained reduction from problem $A$ to problem $B$ with respect to running times $a(n)$ and $b(n)$ implies that an $O(b(n)^{1-\varepsilon})$-time algorithm for $B$ for some $\varepsilon> 0$ can be converted into an $O(a(n)^{1-\delta})$-time algorithm for $A$ for some $\delta> 0$.

Fine-grained complexity took three of the problems from our list above and postulated that the state-of-the-art algorithms cannot be beaten: in the word RAM model of computation with $O(\log n)$-bit words, APSP with polynomially bounded integer weights requires $n^{3-o(1)}$ time (the “APSP hypothesis” [140, 157, 159, 26]), 3SUM on $n$ integers of polynomial size requires $n^{2-o(1)}$ time (the “3SUM hypothesis” [92, 136]), and there is no $\varepsilon> 0$ such that CNF-SAT on $n$ variables and $O(n)$ clauses can be solved in $(2-\varepsilon)^n$ time (the “Strong Exponential Time Hypothesis,” SETH [103, 104, 50]).

Through fine-grained reductions, FGC built a web of relations among problems, and within it, classes of problems that are equivalent to each other (see the surveys [154, 45]). For a large variety of problems for which no improved running times had been obtained in decades, FGC pointed to reasons: they are 3SUM-hard, APSP-hard, SETH-hard, or worse. In other words, to improve the known algorithms, one knew which obstacle needed to be overcome.

For instance, many problems are known to be equivalent to APSP under fine-grained reductions. Either all of the following problems are solvable in truly subcubic time, meaning in $O(n^{3-\varepsilon})$ time for some constant $\varepsilon> 0$ (truly subquadratic is defined similarly), or none of them is: the $(\min,+)$-product $A \star B$ of two $n \times n$ matrices $A$ and $B$ (defined by $(A \star B)[i,j] = \min_k(A[i,k]+B[k,j])$) [87, 19]¹, Negative Triangle and Second Shortest Simple Path [157, 159], Graph Radius [30, 17], Tree Edit Distance [39, 133], and many other problems [157, 159, 138, 30, 17], more of which are listed below, after Figure 1. Perhaps surprisingly, it is also known [66] that solving APSP in truly subcubic time would give a polynomially faster algorithm for directed unweighted APSP as well, and the corresponding hypotheses are in fact equivalent [82]² (albeit conditionally${}^{2}$), even though the longstanding running times of the two problems are different: $\approx n^3$ for APSP and $\approx n^{2.5}$ for directed unweighted APSP (if $\omega= 2$).

Similar reductions and equivalences are also known for 3SUM: all nontrivial 3-linear degeneracy tests are equivalent to 3SUM [73], and many other problems are 3SUM-hard (e.g., many problems in computational geometry [92, 36, 146, 18], and some string problems [6, 26, 116]).

CNF-SAT reduces to Orthogonal Vectors [163], so that unless SETH is refuted, OV requires $n^{2-o(1)}$ time. SETH also implies fine-grained hardness for many exponential-time problems [50].

Both 3SUM and APSP reduce to the Exact Triangle problem (given a tripartite graph with integer edge weights, is there a triangle whose three weights sum to zero?) [158, 157, 159]. Through Exact Triangle, both reduce to online Set Disjointness (from the list above), and even to its offline version, in which the query pairs are given in advance together with the sets [160, 71]; for 3SUM, direct reductions are also known [136, 116, 109, 1]. In graph terms, offline Set Disjointness is the All-Edges Sparse Triangle problem: given a graph with $m$ edges, decide for every edge whether it lies in a triangle.

All-Edges Sparse Triangle thus inherits the hardness of both 3SUM and APSP, and it is in turn the source of the conditional lower bounds for Set Disjointness data structures, dynamic graph problems, and other problems [136, 26, 116, 90]. Many came to regard it as one of the most reliably hard problems in fine-grained complexity, and the evidence seemed good.

All-Edges Sparse Triangle can be solved in $O(m^{3/2})$ time just by listing all triangles [105, 60], and a high-degree/low-degree split combined with matrix multiplication [14] improves this to $O(m^{2\omega/(\omega+1)})$, where $2 \le\omega< 2.372$ is the matrix multiplication exponent.³ Even if $\omega= 2$, this is $m^{4/3}$. On the graphs the reductions produce, which have $n$ vertices of degree about $\sqrt{n}$, this equals the brute-force running time $n^2$ conjectured to be optimal. Matrix multiplication only seemed useful for dealing with dense inputs, whereas these instances ask for a sparse set of outputs, and it seemed plausible that no technique could do better.

However, a reduction from $A$ to $B$, proved in order to show that $B$ is hard, is also an algorithm for $A$ whenever $B$ turns out to be easy. In this paper we give a new algorithm for a lopsided version of All-Edges Sparse Triangle, in which the graph is tripartite and one of the three parts is significantly smaller than the other two. Most known reductions to All-Edges Sparse Triangle can be slightly modified to produce such instances.

As a consequence, all of the problems in the list above (and many more), except for CNF-SAT and Orthogonal Vectors, have polynomially faster algorithms than before:

• 3SUM: deterministic $O(n^{1.9992})$ time for integers of polynomial size (Theorem 22), and for real numbers, a Las Vegas algorithm with $O(n^{1.998})$ expected time whose only operations on real numbers are comparisons, additions, and subtractions (Theorem 35).

• APSP: deterministic $O(n^{2.9995})$ time for polynomially bounded integer weights, and the same bound for the $(\min,+)$-product (Theorem 22), and for real weights, Las Vegas algorithms with $O(n^{2.998})$ expected time, again using only comparisons, additions, and subtractions on real numbers (Theorem 35).

• Directed Unweighted APSP: a deterministic algorithm running in $O(n^{2+\mu-\varepsilon''})$ time for some constant $\varepsilon'' > 0$, where $\mu< 0.5275$ is defined by $\omega(1,\mu,1) = 1 + 2\mu$, and $\omega(1,\mu,1)$ is the exponent of multiplying an $n \times n^\mu$ by an $n^\mu\times n$ matrix. Note that $O(n^{2+\mu}) \le O(n^{2.5275})$ is the longstanding running time of Zwick’s algorithm. This is the first improvement over Zwick’s algorithm that does not come from faster rectangular matrix multiplication, and it refutes the Directed Unweighted APSP hypothesis [66, 68, 84] (Theorem 34).

• Tree Edit Distance: $O(n^{2.9995})$ time for polynomially bounded integer costs via the tight reduction of [39].

• Online Set Disjointness: for $D \le N^{0.12}$, a deterministic data structure whose preprocessing time is polynomially less than $N^2$ and whose query time is polynomially less

than $\sqrt{D}$, and which returns $|S_I \cap T_J|$, so that it solves the counting version as well. For $D \le N^{1/18}$, the preprocessing takes $O(N^2/D^{0.063})$ time and the queries take $O(D^{0.437})$ time (Theorem 3).

The APSP and 3SUM hypotheses are therefore refuted, and so are several other hypotheses of fine-grained complexity. Figure 1 summarizes the problems affected, and in Section 1.2 below we list many more new algorithmic results. All of these results follow from one algorithm, for computing a given sparse set of entries of the product of two thin matrices, which we describe in the next subsection.

We view this as a great success for fine-grained complexity. The hypotheses were made to be tested. The reductions and connections between problems have always been the focus of fine-grained complexity, and these connections persist whether or not the problems turn out to be hard. The web of reductions that was interpreted to imply hardness is exactly what turns a single new algorithm into faster algorithms for a myriad of problems at once. This also opens a number of research directions in fine-grained complexity, algorithm design, and algebraic complexity; we discuss some in Section 6.

All-Edges Sparse Triangle

In this paper we present a simple algebraic technique that solves All-Edges Sparse Triangle in $O(n^{2-\delta})$ time, for a constant $\delta> 0$, on sparse tripartite graphs with two parts of $n$ vertices and a third part of at most $n^{0.12}$ vertices. In matrix language it computes a prescribed sparse set of entries of a dense product of thin matrices. We call a product of an $N \times D$ matrix by a $D \times N$ matrix thin when $D$ is a small power of $N$, and, given a set $W$ of positions of the product, we call these positions, and the entries of the product at them, wanted. Unlike in usual fast matrix multiplication algorithms, the wanted entries are a sparse subset of the full matrix product, and our new algorithm critically exploits this. A matrix multiplication algorithm computes all $N^2$ entries of the product whether they are wanted or not, and thus takes $\Omega(N^2)$ time, whereas ours is faster than this. Our theorem is as follows.

Theorem 1 (see Corollary 26 and Theorem 25). Let $N \ge D^{18}$, let $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$ have entries of absolute value $N^{O(1)}$, and let $W$ be any set of $|W| \le N^2/\sqrt{D}$ positions. The entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in $O(N^2/D^{0.063})$ operations on $O(\log N)$-bit integers. More generally, for every $\varepsilon< 0.1204$ and every $\kappa> 0$ there is a $\gamma> 0$ such that, whenever $D \le N^\varepsilon$ and $|W| \le N^2/D^\kappa$, the task takes $O(N^2/D^\gamma)$ operations.

For comparison, when $|W| = N^2/\sqrt{D}$, computing each wanted entry one by one as an inner product costs $N^2\sqrt{D}$, and computing the whole product with Coppersmith’s rectangular matrix multiplication algorithm [62] costs $N^{2+o(1)}$. Our algorithm is faster than either of these, and uses polynomially less than one operation per entry of the full product.

The algorithm builds on Coppersmith’s algorithm [62] for rectangular matrix multiplication, which makes use of an algebraic identity of Schönhage [143]. The high-level idea is to carefully analyze the inner workings of Coppersmith’s algorithm to see which operations are no longer needed since they are only used to compute output entries outside of $W$. This approach would not work effectively for most known matrix multiplication algorithms,⁴ but it turns out that Schönhage’s identity has sparsity and combinatorial properties we can take advantage of to make this work. See Section 2 for a more detailed overview.

Coppersmith’s paper [62] contains two algorithms. The first algorithm achieves an $n^2 \operatorname{polylog} n$ running time for $n \times n^{0.172} \times n$ matrix multiplication (for matrices with $\operatorname{polylog} n$ bit entries). We do not use this algorithm here, but it has been extensively used before, both in the previous fastest APSP algorithm by Williams [166] and to design the fastest fine-grained algorithms for a number of problems using the polynomial method [164, 30, 29, 66, 16, 8]. Our approach instead builds on the second algorithm in [62], which only achieves an $n^{2+o(1)}$ running time (and not $n^2 \operatorname{polylog} n$) for $n \times n^\alpha\times n$ matrix multiplication for some constant $\alpha> 0$; these $n^{o(1)}$ factors become insignificant in our setting of polynomial savings. Subsequent work on rectangular matrix multiplication has mostly followed this second approach [63, 120, 124, 161] achieving the current bound $\alpha> 0.321$. In order to take advantage of the wanted set, we will modify Coppersmith’s second algorithm quite differently from prior work, in a way which ultimately achieves the worse bound $\alpha> 0.12$.

Computing the wanted entries of a thin matrix product is equivalent to the counting version of Lopsided All-Edges Sparse Triangle (if $X$ and $Y$ are the two biadjacency matrices, then $(XY)[a,b]$ is the number of middle vertices adjacent to both $a$ and $b$), so it solves the detection variant as well. In a tripartite graph with two parts of $n$ vertices and a middle part of $n^\varepsilon$ vertices, our algorithm counts the triangles through each of $|W|$ prescribed pairs of outer vertices in $O(|W| n^{0.437\varepsilon} + n^{2-0.063\varepsilon})$ time for $\varepsilon< 1/18$, which is $O(n^{2-0.063\varepsilon})$ for $|W| \le n^{2-\varepsilon/2}$, and, whenever $|W| \le n^{2-\Omega(1)}$, in truly subquadratic time for every $\varepsilon< 0.1204$ (Corollary 16 and Theorem 25).

Known reductions

It is known that 3SUM [136, 116], Exact Triangle and hence APSP [66], and even their real-valued versions [?], all reduce to Lopsided All-Edges Sparse Triangle. In matrix terms, this detection version is a Boolean matrix problem, and it is also known as offline Set Disjointness [116, 66].

The published reductions are mostly randomized, since the hypotheses were stated for randomized algorithms, but some have been derandomized: the reductions from 3SUM [51, 84], including the reduction of [116] to offline Set Disjointness [84], and the reduction from Exact Triangle to the standard balanced case of All-Edges Sparse Triangle, by Chan and Xu [71]. We show that the reduction from Exact Triangle to the lopsided case can also be made deterministic by hashing the weights modulo a prime chosen with fast matrix multiplication and building the instances as Chan and Xu [71] do (Section 3.2), so we refute both hypotheses deterministically.

With randomization, we refute the hypotheses in their most general form, where the inputs are real numbers. The reduction of Chan, Vassilevska Williams, and Xu [?] from the real-valued problems to sparse triangle problems (which was designed as a hardness argument) uses only additions, subtractions, and comparisons (via Fredman’s trick [87]), so combined with our algorithm it gives Las Vegas algorithms for the real-valued versions of 3SUM and APSP (Section 5.2). Unlike the reductions for integer inputs, which need only the detection version, the reduction of [?] that we use is to the counting version of Lopsided All-Edges Sparse Triangle, which asks for the number of triangles through each prescribed pair, and which our Theorem 5 also solves. Chan, Vassilevska Williams, and Xu [?] also give a reduction to the Boolean version of All-Edges Sparse Triangle, but it is tight only from real APSP, not from real 3SUM or Exact Triangle.

Theorem 2 (Theorems 19, 22, and 35). On a word RAM with $O(\log n)$-bit words, deterministic algorithms solve the following problems, where all numbers in the input are integers of absolute value $n^{O(1)}$:

• Exact Triangle on $n$-vertex graphs in $O(n^{2.9983})$ time,

• APSP on directed $n$-vertex graphs with no negative cycles in $O(n^{2.9995})$ time,

• the $(\min,+)$-product of two $n \times n$ matrices in $O(n^{2.9995})$ time, and

• 3SUM on $n$ numbers in $O(n^{1.9992})$ time.

In the real RAM in which the only operations allowed on real numbers are comparisons, additions, and subtractions, Las Vegas algorithms solve the following problems on real inputs:

• 3SUM on $n$ numbers in $O(n^{1.998})$ expected time,

• APSP in $O(n^{2.998})$ expected time,

• the $(\min,+)$-product in $O(n^{2.998})$ expected time, and

• Exact Triangle in $O(n^{2.998})$ expected time.

Lower bounds are known for the above problems in some restricted models of computation: over the $(\min,+)$ semiring, straight-line programs need $n^3$ operations [111] and path-comparison algorithms need $\Omega(n^3)$ time [112], and 3-linear decision trees for 3SUM need depth $\Omega(n^2)$ [79], [4].

Our algorithms are not in these restricted models, since the known reductions remove the weights (by hashing, or by Fredman’s trick for real inputs), and we then count triangles in unweighted graphs using integer matrix multiplication.

We remark that there already was some indication that some of these models were too restricted, as for instance recent work [113] showed that 6-linear decision trees of depth $O(n\log^2 n)$ suffice to solve 3SUM, and similar strong bounds are known for Exact Triangle. (See the paragraph on unaffected conjectures below for more on this.)

Data structure version, and hinted Online Matrix–Vector. Our algorithm for computing the wanted entries of a thin matrix product can also be converted into a data structure. We strengthen Theorem 1 as follows (Section 4):

Theorem 3 (see Theorem 24 and Corollary 26). Let $N \ge D^{18}$, and let $X \in\mathbb{Z}^{N \times D}$, $Y \in\mathbb{Z}^{D \times N}$ have entries of absolute value $N^{O(1)}$. The pair $(X,Y)$ can be preprocessed deterministically in $O(N^2/D^{0.063})$ time, polynomially less than the size of the product. After this, any single entry $(XY)[I,J]$ can be computed deterministically in $O(D^{0.437})$ time, polynomially less than $\sqrt{D}$. More generally, for every $\varepsilon< 0.1204$ and every $q > 0$ there is a $\gamma> 0$ such that, whenever $D \le N^\varepsilon$, the pair can be preprocessed in $O(N^2/D^\gamma)$ time, after which any single entry can be computed in $O(D^q)$ time.

The bounds of Theorem 1 follow by asking $|W|$ queries to the data structure. This allows us to refute the hinted Online Matrix–Vector (OMv) conjectures of van den Brand, Nanongkai, and Saranurak [155] in the regime of thin hints, i.e., for very rectangular matrices. In the OMv problem, one is given an $n \times n$ Boolean matrix $M$, and then $n$ vectors arrive one at a time, each of which must be multiplied by $M$ before the next one arrives. The OMv conjecture [100] asserts that this requires $n^{3-o(1)}$ total time, and it implies tight lower bounds for many dynamic problems. In the hinted versions, the algorithm gets a hint before the vector arrives. For instance, in $v$-hinted Mv, $M$ is an $n \times t$ matrix and the hint is a $t \times n$ matrix $V$, one of whose columns will be the vector. One can then either compute all of $MV$ by fast matrix multiplication as soon as the hint arrives, or wait and compute one matrix–vector product in $O(nt)$ time, and van den Brand et al. [155] conjecture that no algorithm is polynomially faster than both. Our data structure is faster than both when the hint is thin:

Theorem 4 (Corollary 40). The $v$-hinted Mv, Mv-hinted Mv, and uMv-hinted uMv conjectures of [155] (Conjectures 5.2, 5.7, and 5.12 there) are refuted in the regime of thin hints. For hint dimension $t = n^\tau$ with $0 < \tau< 1/18$, the phase after the hint takes $O(n^{2-0.063\tau})$ time and the phase after the vector $O(n^{1+0.437\tau})$ time, against the conjectured $n^{2-o(1)}$ and $n^{1+\tau-o(1)}$. Correspondingly, the uMv version, whose two hints have dimensions $t_1 = n^{\tau_1}$ and $t_2 = n^{\tau_2}$, fails for $\tau_1 < \tau_2/18$. All three fail, with smaller savings, for every $\tau< 0.1204$ (for the uMv version, $\tau_1 < 0.1204\tau_2$).

The main conditional lower bounds of [155] for dynamic matrix inverse use the conjectures at $\tau\approx0.53$, far from the thin regime, and are unaffected. The bounds that are affected are those for the tail end of the update–query trade-offs with the fastest queries. These had established the conditional optimality of the dynamic matrix inverse algorithms of Sankowski [141] and of [155] when queries are fast, and these no longer have a believable conditional lower bound. Neither do other results whose lower bounds were based on this end of the trade-offs, e.g., for partially dynamic distance oracles [HLS24] and dynamic attention [156] (Section 5.4).

More new algorithms. Figure 1 shows the reductions involved. In addition to the results listed above, composing known reductions with our algorithms for Exact Triangle, 3SUM, APSP, and the $(\min,+)$-product ((2)) gives polynomially faster algorithms for many other problems. We next list some examples. Throughout, the numbers in the input are integers of absolute value polynomial in the input size, $n$ is the number of vertices, matrix rows, or input elements, and our algorithms are deterministic unless stated otherwise.

• The APSP class. Negative Triangle, which asks whether an edge-weighted graph has a triangle of negative total weight, is arguably the simplest problem of the subcubic equivalence class of APSP [157, 159]; the fact that it is equivalent to APSP made many other reductions possible. The APSP subcubic equivalence class also contains Minimum Weight Cycle (with nonnegative weights) [138, 157, 159], Replacement Paths (for each edge of a shortest $s$–$t$ path, the length of a shortest $s$–$t$ path avoiding it) and Second Shortest Simple Path in directed graphs [157, 159], Radius, Median, and (for unique shortest paths) Betweenness Centrality [14], Metricity (does a given matrix satisfy the triangle inequality?) and $(\min,+)$-Product Verification (is $C = A + B$ for given matrices $A$, $B$, and $C$?) [157, 159], Tree Edit Distance [39], Maximum Subarray (the contiguous submatrix of an $n \times n$ matrix with the largest sum of entries) [TT00, BDT16], the Wiener Index (returning the sum of all distances in a graph) [LVW18], computing

the number of witnesses for every entry of the (min, +)-product [68] and a variety of parity versions of many of the mentioned problems [12]. We give the first truly subcubic algorithms for all of these problems, and for Diameter and many other problems which reduce to APSP but may not be known to be equivalent.

• The 3SUM class. The subquadratic equivalence class of 3SUM contains the variant with three input sets $A$, $B$, and $C$ (is there $a \in A$, $b \in B$, $c \in C$ with $a+b=c$?) and GeomBase (given $n$ points with integer coordinates on the three horizontal lines $y=0$, $y=1$, and $y=2$, is there a non-horizontal line through three of them?) [92], All-Numbers 3SUM (for every input number, decide whether it is part of a solution) [157, 159] (the reduction becomes deterministic with the self-reductions of [126, 91], as noted in [84]), Convolution-3SUM (do $x_0,\ldots,x_{n-1}$ satisfy $x_i+x_j=x_{i+j}$ for some $i,j$?) [136, 51], and every nontrivial variant of 3-Linear Degeneracy Testing (are there distinct input numbers $x_1,x_2,x_3$ with $\alpha_1x_1+\alpha_2x_2+\alpha_3x_3=t$, for fixed integers $\alpha_1,\alpha_2,\alpha_3,t$?), such as 3AP, the case $x_1+x_2=2x_3$ (do three of the numbers form an arithmetic progression?) [73]. The counting version #3SUM is also in the class [68, 83]; see the last item below. We give the first truly subquadratic algorithms for all of these problems.

The 3SUM class also contains MonoConvolution, a “colored” version of Boolean convolution: given three integer sequences $a,b,c$ of length $n$, determine for every index $i$ whether there is a $j<i$ with $a[j]=b[i-j]=c[i]$. Lincoln, Polak, and Vassilevska Williams [123], Theorems 7 and 8 show that MonoConvolution is $(n^2,n^{1.5})$ fine-grained equivalent to 3SUM: an $O(n^{2-\varepsilon})$-time algorithm for 3SUM gives an $\widetilde{O}(n^{3/2-\varepsilon/(8-2\varepsilon)})$-time algorithm for MonoConvolution. Our approach gives $O(n^{1.4999})$ time for MonoConvolution, improving on the previous $\widetilde{O}(n^{1.5})$ time of [123].

• The (min, +)-convolution class. The (min, +)-convolution of two sequences $a$ and $b$ of $n$ integers is defined as the sequence $c$ with $c[k]=\min_{i+j=k}(a[i]+b[j])$ for all $k$. This convolution version of APSP is fine-grained reducible to both the (min, +)-product [44] and 3SUM [44, 59]. Problems equivalent to (min, +)-convolution include Superadditivity Testing (is $a[i+j]\ge a[i]+a[j]$ for all $i,j$?), Maximum Consecutive Subsums (for every $k$, the largest sum of $k$ consecutive elements), and Tree Sparsity (a subtree with $k$ vertices of maximum weight in a vertex-weighted tree) [59]. We give the first truly subquadratic algorithms for all of these problems, refuting the (min, +)-convolution hypothesis [44, 59]. (This would have been possible by refuting either the APSP or the 3SUM hypothesis, and we refute both.)

• Knapsack and Approximate Subset Sum. Both 0/1 and Unbounded Knapsack, with $n$ items and capacity $t$, reduce to (min, +)-convolution [59]. So does approximating Subset Sum to within a factor $1-\varepsilon$, which is in fact equivalent to it: an algorithm for (min, +)-convolution in time $T(n)$ gives a randomized approximation scheme in $\widetilde{O}(n+T(1/\varepsilon))$ time [44]. For 0/1 Knapsack, the best known bound in terms of $n$ and the largest item weight $w_{\max}\le t$ is $\widetilde{O}(n+w_{\max}^2)$ [46, 106]. We give the first algorithms running in $\widetilde{O}(n+t^{2-\delta})$ time, and the first approximation scheme running in $\widetilde{O}(n+(1/\varepsilon)^{2-\delta})$ time, for a constant $\delta>0$. Only the algorithm for Unbounded Knapsack is deterministic.

• Zero-, Min-, and Max-Weight $k$-Clique. For a constant $k\ge3$, find $k$ vertices of an edge-weighted graph whose $\binom{k}{2}$ edge weights sum to zero, or to the minimum or maximum possible value. Due to the now classical reduction of Nešetřil and Poljak [132]

$k$-Clique to triangle detection, our new Exact Triangle algorithm immediately implies the first algorithms running in $O(n^{k-\varepsilon})$ time for an $\varepsilon> 0$ (Section 5.3), refuting the weighted $k$-Clique hypotheses [7, 125, 36].

• 3XOR. Given three lists of $n$ vectors in $\mathbb{F}_2^d$, where $d = O(\log n)$, decide whether there are vectors $a$, $b$, and $c$, one from each list, with $a \oplus b \oplus c = 0$, where $\oplus$ denotes the bitwise XOR. This is a well-studied variant of 3SUM [108, 77]. It is not known to reduce to 3SUM, but the same approach works for it as well, so we give the first truly subquadratic algorithm (Remark 23).

• Counting versions. Counting the zero-weight triangles of an Exact Triangle instance, or the solutions of a 3SUM instance, is equivalent to the decision problem [68], and for 3SUM the equivalence also holds under deterministic reductions [83]. Similarly, counting the witnesses of each entry of the $(\min,+)$-product is equivalent to computing the $(\min,+)$-product [68]. We give the first truly subcubic algorithms for #Exact Triangle and for counting the witnesses of the $(\min,+)$-product, and the first truly subquadratic algorithm for #3SUM. See [68, 83] for applications.

For problems that are only known to be 3SUM-hard or APSP-hard, such as many of the geometric problems of [92] and the dynamic problems of [135, 26, 116], we get no faster algorithms, since the reductions go the other way. For instance, it is now open whether one can find three collinear points among $n$ points in the plane in truly subquadratic time.

Relationships to some unaffected conjectures. While we refute various hypotheses, their underlying problems are all similar in nature. We outline how a few other important problems with fine-grained hypotheses differ from 3SUM, Exact Triangle, and APSP, signaling that refuting them might require other techniques.

We begin with CNF-SAT and Orthogonal Vectors (OV). Neither of these problems is known to reduce to All-Edges Sparse Triangle or similar problems. CNF-SAT and Exact Triangle are known to both reduce to certain triangle problems (Matching Triangles and Triangle Collection [28]), but these triangle problems concern the triangles in each of about $n^3$ color triples in dense node-colored graphs, and our techniques do not appear to work there.

Unlike CNF-SAT and Orthogonal Vectors, 3SUM, APSP, and Exact Triangle have long been known to have fast algorithms in certain more powerful models.

3SUM (and more generally $k$-SUM) and also Exact Triangle have very efficient linear decision trees. Following a long line of work [130, 79, 15, 93, 97, 98, 53, 118], it is now known [113] that for every $k \ge3$, $k$-SUM has $2k$-linear decision trees of depth $O(n\log^2 n)$, and Exact Triangle in $m$-edge graphs has 6-linear decision trees of depth $O(m\log^2 m)$. APSP has decision trees of depth $\tilde{O}(n^{2.5})$ [87]. No such results are known for CNF-SAT and OV.

3SUM, Exact Triangle, and APSP also all have fast nondeterministic and co-nondeterministic algorithms. A nondeterministic algorithm for a decision problem verifies a proof for YES, and a co-nondeterministic algorithm verifies a given proof for NO. Carmosino et al. [50] showed that 3SUM can be solved nondeterministically and co-nondeterministically in $\tilde{O}(n^{1.5})$ time, and Exact Triangle and APSP in $\tilde{O}(n^{(3+\omega)/2})$ time. The best known co-nondeterministic algorithms for CNF-SAT and OV are no better than the known deterministic ones, and [50] postulate NSETH, the hypothesis that CNF-SAT requires $2^{n-o(n)}$ time even for co-nondeterministic algorithms. They then prove that if NSETH is true, then there cannot be a deterministic (or zero-error randomized) fine-grained reduction from CNF-SAT or OV to 3SUM, APSP, or Exact Triangle. In fact, our algorithm shows that if there were a deterministic fine-grained reduction from CNF-SAT or OV to 3SUM, APSP, or Exact Triangle, then SETH (and not merely NSETH) is false.

Interestingly, the issues posed by co-nondeterminism disappear in the Merlin–Arthur world where randomization is allowed. CNF-SAT [165], 3SUM, and APSP [3] all have vastly improved Merlin–Arthur protocols: the verification time for CNF-SAT (indeed for counting its satisfying assignments) is $2^{n/2}\operatorname{poly}(n)$, for 3SUM it is $\widetilde{O}(n)$, and for APSP it is $\widetilde{O}(n^2)$.

While we obtain a polynomial improvement for 3SUM, we are unable to do so for $k$-SUM for $k \ge4$. 3SUM reduces to $k$-SUM, but no reduction in the other direction is known, and although $k$-SUM reduces to exact-weight subgraph problems [21, 23], whose $n^k$ baseline we improve (Section 5.3), those reductions can only give $k$-SUM algorithms far slower than the $n^{\lceil k/2\rceil}$ baseline.

The same holds for $k$-XOR for $k \ge4$ and for the 3SUM-Indexing conjecture [115], [115], [88] . Ultimately, we believe that one good reason for the lack of progress on $k \ge4$ is that $k$-SUM for $k \ge4$ is not known to have fine-grained self-reductions. This has been a major obstacle in designing reductions from $k$-SUM and is also a major obstacle in designing algorithms for it. (The same goes for $k$-XOR for $k \ge4$.)

Combinatorial algorithms. The new algorithms are algebraic and potentially impractical in their current form: the constants hidden in the $O(\cdot)$ are enormous, and the exponents can likely be improved. One might ask whether there are practical or “combinatorial” algorithms for these problems, similar to how one talks about combinatorial Boolean Matrix Multiplication (BMM). The known algorithms for combinatorial BMM, including the Four Russians [9], the $n^3/\operatorname{polylog} n$ bounds of [34], [54], [169], and the recent $n^{3/2}/2^{\Omega((\log n)^{1/7})}$ bound of Abboud, Fischer, Kelley, Lovett, and Meka [11], save only a subpolynomial factor, and no truly subcubic combinatorial algorithm is known. We do not know whether such algorithms exist. Most FGC reductions are combinatorial, so one could refocus the hypotheses to be about combinatorial algorithms [157, 159, 26].

Other parameter regimes. Even without considering combinatorial algorithms, one might also ask whether All-Edges Sparse Triangle can be solved faster in every regime of how large the three parts are. Our technique (due to the identity we use) needs the middle part to have at most $n^{0.12}$ vertices, whereas prior reductions have focused on different regimes.

For instance, the balanced sparse case of All-Edges Sparse Triangle introduced by Pătraşcu [135] has as input a tripartite graph with roughly equal parts and $m$ edges overall. Most known fine-grained reductions produce such balanced instances and give conditional lower bounds of $m^{4/3-o(1)}$.

Another case of interest is where the tripartite graph has two parts of size $n$ and one part of size $\sqrt{n}$ and $O(n^{3/2})$ edges overall. Solving this version in $O(n^{2-\varepsilon})$ time, even only for triangle detection (and not even for the all-edge version), would have consequences for girth approximation in undirected unweighted graphs [139].

Even faster algorithms. Our focus in this paper is to give the simplest presentation of our new algorithm for thin matrix products, and show how known reductions to All-Edges Sparse Triangle combine with this to give all the polynomial speedups mentioned above. In particular, over the course of this project, with Claude we generated some algorithms whose exponents are slightly better than the ones we present here, but which are significantly more complicated. We intentionally did not include them in this paper. See also the discussion in Section 6 and footnote 10.

Organization. In Section 2 we give the main matrix algorithm. In Section 3 we derive the algorithms for Exact Triangle, 3SUM, and APSP from it, by using or slightly modifying known reductions. In Section 4 we give the data structure version of our matrix algorithm, proving Theorems 1 and 3. In Section 5 we present reductions to achieve the other algorithmic results discussed above. Finally, in Section 6 we discuss future directions.

Quickly computing certain entries of a thin matrix product

We now present our main new algorithm, for computing a subset of the output entries of a thin matrix product. We call such a matrix product thin since the input matrices are highly rectangular. We work in the standard word RAM model with $O(\log N)$-bit words, and thus count operations on $O(\log N)$-bit integers.

Theorem 5. Let $D \ge4$ be a power of four and $N \ge D^{18}$. Given as input matrices $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$, whose entries are integers of absolute value at most $N^{O(1)}$, as well as a set $W$ of at most $N^2/\sqrt{D}$ positions of an $N \times N$ matrix, the entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in time $O(N^2\log^2 D/D^{1/18})$.

We call the positions in $W$, and the entries of $XY$ at them, wanted. We emphasize that $X$, $Y$, and $XY$ may be arbitrary dense matrices, and no structure is assumed about the set $W$ other than its size $|W| \le N^2/\sqrt{D}$.

We have aimed to present this section in a way that is accessible to nonexperts. Some details in Theorem 5, such as the exponent 18, are chosen to simplify the presentation in this section while still proving a version sufficient to refute the 3SUM and APSP hypotheses. Later, in Section 4, we present a stronger data structure version that has improved parameter trade-offs.

Notably, while the language of tensors and bilinear complexity naturally describes much of our algorithm and the prior work it builds on, we deliberately avoid that language here. This is in part so that a reader unfamiliar with that background may still follow, but also because our arguments here will focus on sparsity and combinatorial properties of our algorithms and algebraic identities that are often intentionally abstracted away in tensor language.

Sections 2.1–2.3 present an exposition of the relevant prior work that we directly build on, especially Strassen [147], Schönhage [143], and Coppersmith [62], and no prior understanding of those works is assumed. The main new idea of this paper is then presented in Section 2.4.

Technique overview. Strassen’s algorithm multiplies matrices by applying a small identity (seven products for multiplying $2 \times2$ matrices) recursively. In that algorithm, the multiplications happen at the leaves of a recursion tree, and the additions that prepare the leaves (the “encoding phase”) and assemble the result (the “decoding phase”) are comparatively cheap. We do the same with an identity of Schönhage that computes, simultaneously, two small rectangular matrix products. Applied at $L$ levels of recursion, the identity computes $2^L$ matrix products at once. In order to multiply $X$ and $Y$ using this, we cut them into smaller blocks, and then simultaneously multiply many pairs of blocks using Schönhage’s identity (Section 2.3).

The new idea (Section 2.4) is to visit only the leaves of the recursion tree that the wanted entries from $W$ need. A priori, it is not clear that this reduces the number of leaves by much, or that the encoding and decoding phases stay cheap. Both nonetheless turn out to hold, for two reasons.

First, we observe that the additions in the encoding phase that prepare the leaves depend on $X$ alone or on $Y$ alone (as is the case in all tensor rank decompositions). Since each $X$ block is multiplied with every $Y$ block, we only need to compute these additions once for each $X$ block and may share them among many calls to our algorithm (and similarly for $Y$). One reason for the requirement $N \ge D^{18}$ is to ensure that this total cost will be negligible (Section 2.4.1).

Second, we prune the recursion tree, and make a recursive call only when one of its outputs is wanted by $W$. This way, the running time is governed by the union of the sets of leaves that the wanted entries need. Most of the analysis goes into bounding this union (Section 2.4.3). We will see that every output entry has a private leaf that no other entry needs, but that all its other leaves are shared by many entries. We will characterize each leaf by how often one needs to pick each of the two matrix products from Schönhage’s identity to get to it in the recursion tree, and find a trade-off: choosing one of the matrix products leads to a relatively small number of leaves, while choosing the other increases how many output entries share each leaf. Balancing the two, we will find that a set of $\frac{N^2}{\sqrt{D}}$ wanted entries only needs about $\frac{N^2}{D^{1/18}}$ leaves in all, and paying $O(\log^2 D)$ operations per leaf gives Theorem 5.

Strassen’s recursive algorithm

We use the recursive approach for multiplying matrices introduced by Strassen [147]. We begin by recalling his algorithm for intuition. The key behind his algorithm is an identity that shows how to multiply two $2 \times2$ matrices $A$ and $B$ using only seven multiplications instead of eight. It forms seven products, each of a linear combination of the entries of $A$ with a linear combination of the entries of $B$, and then obtains each entry of $AB$ as a linear combination of those seven products. Recall that, writing $A = (a_{ij})$ and $B = (b_{ij})$, the seven products are

and the product is

To use this identity to multiply $N \times N$ matrices $A,B$ for larger $N$, we interpret $A$ and $B$ as $2 \times2$ block matrices (whose entries are $N/2 \times N/2$ matrices), then multiply them as in the identity, i.e., we form seven pairs of linear combinations of the blocks of $A$ and of $B$, recursively multiply each pair, and then assemble the four blocks of $AB$ from the seven results. When $N = 2^L$, the recursion tree (Figure 2) has $L$ levels, and $7^L$ leaves, one for each choice of one of the seven products at each level. The multiplications take place at the leaves, meaning there are $7^L = N^{\log_2 7} < N^{2.81}$ multiplications in total. We call the additions on the way down, which prepare the leaves and are computed from $A$ alone and from $B$ alone, the encoding, and those on the way up, which assemble the product, the decoding.

To calculate the total cost of Strassen’s algorithm, we also need to count the additions and subtractions at the internal vertices of the recursion tree. At depth $k$ of the recursion tree, there are $7^k$ internal vertices. Both on the way down and on the way up the tree, the algorithm must compute linear combinations of $2^{L-k} \times2^{L-k}$ matrices at each of these vertices, which uses $O(4^{L-k})$ operations. The total across the whole tree is therefore $\sum_{k=0}^{L-1} 7^k4^{L-k} = O(7^L)$, a geometric series dominated by its last term, so the additions only increase the running time by a constant factor. This phenomenon is prevalent in recursive matrix multiplication algorithms, although in the new algorithm we develop in Section 2.4 below, this same argument will no longer apply, and we will need to count the additions and subtractions more carefully.

Schönhage’s identity for an inner product and an outer product

While we will still use Strassen’s recursive approach here, we will replace his identity with another one due to Schönhage [143]. In this subsection we give the notation to introduce this identity and explain its usefulness.

It will be convenient to write such identities as trilinear polynomials.⁵ As an example, we first show how to write Strassen’s identity from Section 2.1 in this way. Introduce four output variables $c_{11}, c_{12}, c_{21}, c_{22}$, one for each entry of the product $AB$, and multiply each of the seven products $S_1, \ldots, S_7$ of Strassen’s identity above by the signed sum of the output variables of the entries where it appeared in the decoding phase. For instance, $S_1$ appears in the formulas for $(AB)_{11}$ and $(AB)_{22}$, so it is multiplied by $c_{11} + c_{22}$, and $S_5$ appears in the formula for $(AB)_{11}$ with a minus sign and in the formula for $(AB)_{12}$, so it is multiplied by $c_{12} - c_{11}$. Strassen’s identity then reads

where $S_1, \ldots, S_7$ are the products of two linear forms from (1). Now the coefficient of $c_{ik}$ on the right-hand side is the entry $(AB)_{ik} = a_{i1}b_{1k} + a_{i2}b_{2k}$ that we want, and its coefficient on the left-hand side is the combination of the seven products that computes it, such as $S_1 + S_4 - S_5 + S_7 for $c_{11}$. This polynomial identity thus records the linear combinations and products to be formed, as well as the combinations that assemble the outputs from them.

We will next present Schönhage’s identity in the same way. Recall that the outer product of two vectors $(x_1,x_2,x_3)$ and $(y_1,y_2,y_3)$ consists of the nine numbers $x_i y_j$, and the inner product of two vectors $(p_{11},p_{12},p_{21},p_{22})$ and $(q_{11},q_{12},q_{21},q_{22})$ is the single number $\sum_{i,j=1}^{2}p_{ij}q_{ij}$. These are special cases of matrix products: the outer product is multiplying a $3 \times1$ matrix times a $1 \times3$ matrix, and the inner product is multiplying a $1 \times4$ matrix times a $4 \times1$ matrix. (We index the $p$ and $q$ variables with $\{1,2\} \times\{1,2\}$ rather than $\{1,2,3,4\}$ to set up nicer notation for $\hat{p}$ and $\hat{q}$ in Lemma 6 below.)

Separately computing the inner product and outer product would use $9+4=13$ multiplications. However, Schönhage [143], §6, eq. (6.3), with $k=n=3$ showed that they can be computed (along with some error terms that will turn out to be harmless) using only 10 multiplications. The fact that this is possible even when the inner and outer products are on disjoint sets of variables is surprising, and refuted a border rank version of the “direct sum conjecture for tensor rank” [148].

We call $x_1,x_2,x_3,p_{11},p_{12},p_{21},p_{22}$ the seven left variables and $y_1,y_2,y_3,q_{11},q_{12},q_{21},q_{22}$ the seven right variables. We will also use ten output variables $z_{ij}$ (for $1 \leq i,j \leq3$) and $z_0$, one for each of the ten quantities we aim to compute, namely, the coefficients of the output variables in the polynomial

As a notational helper to state the identity, we introduce two $3 \times3$ matrices $\hat{p}$ and $\hat{q}$ of linear forms. These are formed by starting with the $2 \times2$ matrices of the $p$’s and of the $q$’s, then extending them to $3 \times3$ matrices so that every column of $\hat{p}$ and every row of $\hat{q}$ sums to zero:

Notice that $\sum_{i,j=1}^{3}\hat{p}_{ij}\hat{q}_{ij}=\sum_{i,j=1}^{2}p_{ij}q_{ij}$ is the desired inner product.

Schönhage’s identity has ten terms, which we name $P_{ij}$ (for $1 \leq i,j \leq3$, corresponding to the entries of $\hat{p}$ and $\hat{q}$) and $P_0$. Each term $\lambda$ is the product of a linear form $\varphi_\lambda$ in the left variables, a linear form $\psi_\lambda$ in the right variables, and a linear form $\chi_\lambda$ in the output variables. These linear forms are as follows:

Lemma 6. Summing over the ten terms,

Proof. We compute the coefficient of each output variable on the left-hand side. The coefficient of $z_{ij}$ comes only from the term $P_{ij}$, and it is $\varphi_{P_{ij}}\psi_{P_{ij}}=(x_i+\hat{p}_{ij})(y_j+\hat{q}_{ij})=x_i y_j+(x_i\hat{q}_{ij}+\hat{p}_{ij}y_j+\hat{p}_{ij}\hat{q}_{ij})$. The coefficient of $z_0$ is $\sum_{i,j=1}^{3}(x_i+\hat{p}_{ij})(y_j+\hat{q}_{ij})-(x_1+x_2+x_3)(y_1+y_2+y_3)$. In this expression, the products $x_i y_i$ cancel, and the cross terms sum to $\sum_{i=1}^{3} x_i \sum_{j=1}^{3} \hat{q}_{ij} + \sum_{j=1}^{3} y_j \sum_{i=1}^{3} \hat{p}_{ij} = 0$ since the rows of $\hat{q}$ and the columns of $\hat{p}$ sum to zero. What remains is $\sum_{i,j=1}^{3} \hat{p}_{ij}\hat{q}_{ij} = \sum_{i,j=1}^{2} p_{ij}q_{ij}$, as desired.

Lemma 6 is the key identity we will use. It computes our desired $G$ plus an error polynomial $E$. Before moving on to our algorithm, we make two structural observations about the identity that we will use below.

Our first observation will be used to handle the error $E$. We call the variables of the inner product, $p$, $q$, and $z_0$, inner, and those of the outer product, $x$, $y$, and $z_{ij}$, outer. Notice that every monomial of $G$ consists only of inner variables ($p_{ij}q_{ij}z_0$) or only of outer variables ($x_i y_j z_{ij}$). Meanwhile, every monomial of $E$ has an outer output variable together with an inner left or right variable. We will use this observation to show that $E$ is harmless in our algorithm in Lemma 9 below.⁶

Our second observation concerns the pattern of output variables in our identity, and will be used to take advantage of the sparsity of the set $W$ of wanted entries in our final algorithm. We say that a term $\lambda$ contributes to an output variable $z$ if $z$ appears in $\chi_\lambda$. We see that every term contributes to $z_0$, and only $P_{ij}$ contributes to $z_{ij}$. (The decoding layer of Figure 3 below will also illustrate this.) The fact that only one term contributes to each $z_{ij}$ will eventually help us to reduce how many products we need to compute when we only need certain output entries. Since the $\chi_\lambda$ forms are all sums of output variables (without any minus signs or other coefficients), it follows that for each output variable $z$, its coefficient on the left-hand side of the identity is $\sum_\lambda\varphi_\lambda\psi_\lambda$, where the sum is over the terms $\lambda$ contributing to $z$.

Rectangular matrix multiplication using Schönhage's identity

In this subsection, we apply Schönhage’s identity recursively to derive an algorithm for rectangular matrix multiplication (computing all the output entries, not just a small subset). Our approach is similar to Coppersmith’s algorithm [62], although we will make modifications that worsen some parameters of the algorithm in exchange for preserving some of the combinatorial and structural properties of Schönhage’s identity that we will use below.

The recursion

Schönhage’s identity computes two matrix products, the outer product and the inner product, on disjoint sets of variables. We will see that applying it to $L$ levels of recursion therefore involves $2^L$ matrix products, one for each way of choosing the outer or the inner product at each level. A choice with the inner product at $m$ of the levels, and the outer product at the other $L-m$ levels, corresponds to the product of a $3^L \times4^m$ matrix and a $4^m \times3^{L-m}$ matrix. (We will focus on a setting of $m$ much smaller than $L$ so that these are thin matrix products.) In this subsection, we set up notation for the arrays that the recursion works on and state the recursion, then in Section 2.3.2 we will describe precisely what it computes.

The inputs and outputs of our algorithm will be arrays that contain the entries of many different matrices concatenated together. (Later in Section 2.3.4 we will use this algorithm as a subroutine to multiply a single pair of matrices.) We will index the entries of these arrays by strings of $L$ variables of Schönhage’s identity. Whether each variable is an inner or outer variable will tell us which of the $2^L$ matrix products the array entry corresponds to, and the specific variables will tell us which entry of the matrix it corresponds to.

A left string of length $L$ is a string $u = u_1u_2\cdots u_L$ of $L$ left variables, and we call $u_\ell$ its variable at level $\ell$. An array $a$ on the left strings of length $L$ assigns an integer $a[u]$ to each $u$ that is a left string of length $L$. The array $a$ has $7^L$ entries in total (since there are 7 left variables). We define right strings, output strings, and arrays on them in the same way. For instance, an array on the output strings of length $L$ has $10^L$ entries (since there are 10 output variables).

For a term $\lambda$ and a left variable $s$, let $\varphi_\lambda(s) \in\{0,\pm1\}$ denote the coefficient of $s$ in $\varphi_\lambda$, and define $\psi_\lambda(t)$ similarly for right variables $t$. Recall from Schönhage’s identity that, for each term $\lambda$, at most three of the seven values $\varphi_\lambda(s)$ are nonzero since each $\varphi_\lambda$ has at most three monomials, and similarly for the seven values $\psi_\lambda(t)$.

Just as Strassen’s algorithm cuts matrices into blocks, we cut an array $a$ on strings of length $L \ge1$ into slices according to the variable at level 1: the slice of $a$ at a left variable $s$ is the array $a_s$ on the strings of length $L-1$ given by $a_s[u'] := a[su']$. Thus, $a$ consists of its seven slices, just as a block matrix consists of its blocks. We similarly define the slices $b_t$ of an array $b$ on right strings, and $c_z$ of an array $c$ on output strings.

The recursion Full$(a,b)$, for arrays $a,b$ on the left and right strings of length $L$, returns an array on the output strings of length $L$:

(1) If $L = 0$, return the number $ab$.

(2) For each term $\lambda$ of Schönhage’s identity, form $A_\lambda:= \sum_s \varphi_\lambda(s)\,a_s$ and $B_\lambda:= \sum_t \psi_\lambda(t)\,b_t$.

(3) For each term $\lambda$, recursively compute $C_\lambda:= \operatorname{Full}(A_\lambda,B_\lambda)$.

(4) Return the array $c$ whose slices are $c_{z_{ij}} := C_{\rho_{ij}}$ and $c_{z_0} := \sum_\lambda C_\lambda$.

For $L = 1$, Full is one application of Schönhage’s identity: step (2) evaluates the linear forms $\varphi_\lambda$ and $\psi_\lambda$ at the seven numbers $a[s]$ and the seven numbers $b[t]$, step (3) multiplies the two values for each $\lambda$, and step (4) sums the products into the ten outputs as the forms $\chi_\lambda$ prescribe. In the terms of Section 2.1, step (2) is the encoding and step (4) the decoding.

The calls of Full form a recursion tree. We refer to the call reached by choosing the terms $\tau_1,\ldots,\tau_k$ at levels $1,\ldots,k$ as the vertex $\tau_1\cdots\tau_k$, a string of $k$ terms. It is at depth $k$, and its input arrays are on strings of length $L-k$. The root is the vertex with $k = 0$, and the leaves, the vertices $\tau= \tau_1\cdots\tau_L$ at depth $L$, are the $10^L$ multiplications performed in step (1).

Let us now count the operations performed by this algorithm. At depth $k$ there are $10^k$ vertices (since Schönhage’s identity has 10 terms), each with two input arrays of $7^{L-k}$ entries and an output array of $10^{L-k}$ entries. Step (2) takes $O(7^{L-k})$ operations at a vertex, and step (4) takes $O(10^{L-k})$ operations. Hence, Full uses $\sum_{k=0}^{L} 10^k \cdot O(7^{L-k}+10^{L-k}) = O(L\cdot10^L)$ operations in total, which is not much more than the $10^L$ multiplications. Of this, the encodings, step (2), account for only $\sum_{k=0}^{L} 10^k \cdot O(7^{L-k}) = O(10^L)$ operations, and the decodings, step (4), for the rest.

Unraveling the recursion computation

In this subsection, we carefully work out which two numbers are multiplied at the leaf $\tau$, and what is returned by $\operatorname{Full}(a,b)$. For a leaf $\tau$ let

where the sums are over all left strings $u$ and all right strings $v$ of length $L$. We extend the relation “contributes to” of Section 2.2 from terms to vertices: a vertex $\tau_1\cdots\tau_k$ contributes to an output string $w$ if $\tau_\ell$ contributes to $w_\ell$ at every level $\ell\le k$. Thus the leaf $\tau$ contributes to $w$ if $\tau_\ell=P_{ij}$ at every level $\ell$ where $w_\ell=z_{ij}$ (but $\tau_\ell$ may be arbitrary at the levels where $w_\ell=z_0$). Let $\operatorname{Mult}(a,b)$ be the array on the output strings of length $L$ given by

Lemma 7. At the leaf $\tau$, $\operatorname{Full}(a,b)$ multiplies $\Phi_\tau(a)$ by $\Psi_\tau(b)$, and it returns $\operatorname{Mult}(a,b)$.

Proof. By induction on $L$ (see Figure 3). For $L=0$, the only leaf is the empty string $\tau$, and $\Phi_\tau(a)=a$ and $\Psi_\tau(b)=b$, as the products over the levels are empty, so step (1) multiplies them and returns $ab=\operatorname{Mult}(a,b)$. For $L\ge1$, step (2) is the encoding: since $a_s[u']=a[su']$, for every leaf $\lambda\tau'$ we have

and similarly $\Psi_{\lambda\tau'}(b)=\Psi_{\tau'}(B_\lambda)$. Thus, by induction, the leaf $\lambda\tau'$ multiplies $\Phi_{\lambda\tau'}(a)$ by $\Psi_{\lambda\tau'}(b)$, and $C_\lambda[w'] = \sum_{\tau'} \Phi_{\lambda\tau'}(a)\Psi_{\lambda\tau'}(b)$, summed over the $\tau'$ contributing to $w'$. Step (4) is the decoding: the leaves contributing to $zw'$ are the $\lambda\tau'$ with $\lambda$ contributing to $z$ and $\tau'$ contributing to $w'$, so $c[zw'] = \sum_{\lambda\text{ contributing to }z} C_\lambda[w']$ agrees with (2).

Since one application of the identity computes the coefficients of $G+E$, we expect that $L$ recursive applications should compute products of $L$ such coefficients. For a left variable $s$, a right variable $t$, and an output variable $z$, let

the coefficient of the monomial $stz$ on the left-hand side of the identity. By Lemma 6, it is also the coefficient of $stz$ in $G+E$, and we know that $\gamma(s,t,z)\in\{0,\pm1\}$.

Lemma 8. For all arrays $a,b$ and every output string $w$,

the sum over all left strings $u$ and right strings $v$.

Proof. By the definitions of $\Phi_\tau$ and $\Psi_\tau$, the coefficient of $a[u]b[v]$ in (2) is $\sum_\tau\prod_{\ell=1}^{L}\varphi_{\tau_\ell}(u_\ell)\psi_{\tau_\ell}(v_\ell)$, summed over the leaves $\tau$ contributing to $w$. These are the leaves whose term at level $\ell$ contributes to $w_\ell$, independently for each $\ell$, so this sum of products is the product of sums

by (3).

We will use both descriptions of $\operatorname{Mult}$: Lemma 8 to read matrix products off it (Section 2.3.3), and the sum over leaves in (2) to skip the leaves that do not contribute to the wanted entries (Section 2.4).

Batch computation of multiple matrix products

We now explain the claim from the beginning of Section 2.3.1, that our recursive algorithm is simultaneously computing many matrix products. We will also show that the error $E$ is harmless, using the fact that each of its monomials uses an outer output variable, but at least one inner input variable.

The inner set of a left string is the set of levels at which it has a $p$ (as opposed to an $x$). It will tell us which of the $2^L$ matrix products the string belongs to. Similarly, the inner set of a right string is the set of levels at which it has a $q$ (as opposed to a $y$), and the inner set of an output string is the set of levels at which it has $z_0$ (as opposed to a $z_{ij}$). The inner part of a string is the string of its variables at the levels of its inner set, in the order of the levels, and its outer part is the string of its variables at the other levels. A string is determined by its inner set, its outer part, and its inner part.

Fix $m\ge1$, and let $N_0:=3^{L-m}$ and $D:=4^m$. Of the $2^L$ matrix products, one for each set $Q\subseteq\{1,\ldots,L\}$ of levels at which the inner product is chosen, we will use only the $\binom{L}{m}$ products with $|Q| = m$, and set the inputs of the others to zero. As we now make precise, the product for such a $Q$ multiplies an $N_0 \times D$ matrix by a $D \times N_0$ matrix: the $L-m$ outer levels, with three left and three right outer variables each, give its $N_0 = 3^{L-m}$ rows and $N_0$ columns, and the $m$ inner levels, with four left and four right inner variables each, give the $D = 4^m$ terms of the sum in each of its entries. (To prove Theorem 5, we will ultimately set $L = 19m$, a choice we will explain in Section 2.4.3 below.)

We will index the rows and columns of the matrices in these products by strings of variables. For an $N_0 \times D$ matrix, we index its $N_0$ rows by the $N_0$ strings of $L-m$ outer left variables, and its $D$ columns by the $D$ strings of $m$ inner left variables. Similarly, for a $D \times N_0$ matrix, we index its $D$ rows by the $D$ strings of $m$ inner right variables, and its $N_0$ columns by the strings of $L-m$ outer right variables. To multiply an $N_0 \times D$ matrix $A$ by a $D \times N_0$ matrix $B$, we match the column of $A$ indexed by the string $\pi= p_{i_1j_1}p_{i_2j_2}\cdots p_{i_mj_m}$ with the row of $B$ indexed by the string $\pi' = q_{i_1j_1}q_{i_2j_2}\cdots q_{i_mj_m}$ with the same indices, so that the entry of $AB$ in row $r$ and column $c$ is

where the sum is over all the $D$ strings $\pi$ of $m$ inner left variables. Thus, the rows of the product are indexed by the strings of outer left variables and its columns by the strings of outer right variables.

Let $Q$ be a subset of $\{1,\ldots,L\}$ of size $m$. As above, the outer part of a left string $u$ with inner set $Q$ is a row of an $N_0 \times D$ matrix and its inner part is a column, so the left strings with inner set $Q$ index the entries of such a matrix. For an $N_0 \times D$ matrix $X_Q$ we write $X_Q[u]$ for its entry at the row and the column of $u$. Similarly, the right strings $v$ with inner set $Q$ index the entries of a $D \times N_0$ matrix $Y_Q$, the row given by the inner part of $v$ and the column by the outer part, and we write $Y_Q[v]$ for this entry. Finally, the output strings with inner set $Q$ index the entries of the product $X_QY_Q$: if the variables of $w$ at the levels outside $Q$ are $z_{i_1j_1},\ldots,z_{i_{L-m}j_{L-m}}$, in the order of the levels, then the row of $w$ is the string $x_{i_1}\cdots x_{i_{L-m}}$ and its column is the string $y_{j_1}\cdots y_{j_{L-m}}$, and we write $(X_QY_Q)[w]$ for the entry there (Figure 4).

We can now check that Full, which returns Mult$(a,b)$, really computes these products. For each of the $K := \binom{L}{m}$ subsets $Q$ of $\{1,\ldots,L\}$ of size $m$, let $X_Q$ be any $N_0 \times D$ matrix and $Y_Q$ be any $D \times N_0$ matrix. Let the input arrays $a$ and $b$ be these matrices laid side by side, by setting

for all the left strings $u$ and the right strings $v$ whose inner set $Q$ has exactly $m$ elements, and $a[u] := 0$, $b[v] := 0$ at every other string.

Lemma 9. For the input arrays $a,b$ above and every output string $w$ whose inner set $Q$ has exactly $m$ elements,

Proof. By Lemma 8, $\operatorname{Mult}(a,b)[w]$ is the sum, over all left strings $u$ and right strings $v$, of $\gamma(u_1,v_1,w_1)\cdots\gamma(u_L,v_L,w_L)a[u]b[v]$, where $\gamma(s,t,z)$ is the coefficient of $stz$ in $G+E$. We determine which summands can be nonzero. First, $a[u]b[v]\ne0$ requires the inner sets of $u$ and $v$ to have exactly $m$ elements. Next, at a level $\ell\in Q$ we have $w_\ell=z_0$, and the only monomials of $G+E$ that contain $z_0$ are $p_{ij}q_{ij}z_0$, with coefficient $1$. So $\gamma(u_\ell,v_\ell)=1$ for some $i,j$. Thus $u$ and $v$ have an inner variable at every level of $Q$, and their inner sets are therefore exactly $Q$. Finally, at a level $\ell\notin Q$ we have $w_\ell=z_{ij}$ for some $i,j$, and $u_\ell$ and $v_\ell$ are outer. The left string $u$: $x_1$ $p_{11}$ $x_2$ $x_3$ $p_{22}$ $x_1$

outer part $x_1x_2x_3x_1$: the row of $u$ in $X_Q$, inner part $p_{11}p_{22}$: its column, so $u$ indexes an entry of $X_Q$

right string $v$: $y_3$ $q_{11}$ $y_1$ $y_3$ $q_{22}$ $y_2$

inner part $q_{11}q_{22}$: the row of $v$ in $Y_Q$, outer part $y_3y_1y_3y_2$: its column, so $v$ indexes an entry of $Y_Q$

output string $w$: $z_{13}$ $z_0$ $z_{21}$ $z_{33}$ $z_0$ $z_{12}$

row $x_1x_2x_3x_1$ and column $y_3y_1y_3y_2$, read off from its $z_{ij}$, so $w$ indexes an entry of $X_QY_Q$

levels $\ell= 1,\ldots,6$ (shaded: the inner set $Q = \{2,5\}$)

only monomial of $G + E$ that contains $z_{ij}$ and no inner variable is $x_iy_jz_{ij}$, with coefficient 1, so $(u_\ell,v_\ell) = (x_i,y_j)$.

Hence, the summands that remain are those where $u$ has inner set $Q$, the row $r$ of $w$ as its outer part, and some inner part $\pi$, and where $v$ has inner set $Q$, the column $c$ of $w$ as its outer part, and the inner part $\pi'$ with the same indices as $\pi$. There is one such summand for each of the $D$ strings $\pi$, and its coefficient is 1. Their sum is $\sum_\pi X_Q[r,\pi]Y_Q[\pi',c] = (X_QY_Q)[w]$.

To summarize: one run of Full, which uses $10^L$ multiplications, computes all $K$ of the products $X_QY_Q$. The $K$ products can be completely independent: we are free to plug in any choices for $X_Q$ and $Y_Q$ and the algorithm will give us all $K$ products. In the next subsection we will end up plugging in the same matrix for different copies $X_Q$ (and $Y_Q$).

Tiling the $N \times D \times N$ product by products of shape $N_0 \times D \times N_0$

We now use our algorithm, which performs $K$ simultaneous multiplications of an $N_0 \times D$ matrix with a $D \times N_0$ matrix, in order to do a single multiplication of an $N \times D$ matrix $X$ with a $D \times N$ matrix $Y$. For our parameter setting, $N \gg N_0$, so we can do this with a simple tiling approach. Tiling enlarges only the outer dimensions, so it handles a smaller inner dimension than Coppersmith’s algorithm from the same identity (at best $D \approx N^{0.1204}$ rather than $N^{0.1402}$ [62]). We use it because each entry of $XY$ is then a single output of a single run of Full, rather than a linear combination of many outputs as in Coppersmith’s algorithm. We will take advantage of this in our algorithm below (Section 2.4).

We cut $X$ into row blocks of $N_0$ consecutive rows and $Y$ into column blocks of $N_0$ consecutive columns. Thus, $XY$ consists of all the products of a row block by a column block, each an $N_0 \times D$ by $D \times N_0$ product, and one run of Full computes $K$ such products at once. We use each run on a $K_0 \times K_0$ grid of block products, where $K_0 := \lfloor\sqrt{K}\rfloor$: we group the row blocks into bands of $K_0$ consecutive blocks, and similarly the column blocks, and we call a row band together with a column band a tile (Figure 5). Each tile is one run of Full. We fix $K_0^2 \le K$ distinct subsets of $\{1,\ldots,L\}$ of size $m$, one for each block product of the grid, the same in every tile. In a tile, let $X_Q$ and $Y_Q$ be the row block and the column block of the block product with subset $Q$, and let $X_Q = Y_Q = 0$ for the other subsets $Q$. Then one run of Full computes all $K_0^2$ block products of the tile (Lemma 9), and running it on every tile computes all of $XY$.

We pair the blocks in a grid like this, rather than arbitrarily, so that the left input array of a tile depends only on its row band and the right one only on its column band. This will let us share their encodings among tiles in Section 2.4.1 below. (To make the bands fit, we first pad $N$ to a multiple of $K_0N_0$ with zero rows of $X$ and zero columns of $Y$, which at most doubles $N$ since $N \ge K_0N_0$, as we will see in Section 2.4.4.)

Before moving on, let us compute the number of operations performed by this algorithm. Let $M := KN_0^2$ denote the number of output entries among all $K$ products of a tile. Since $K_0 \ge\sqrt{K}/2$, there are $(N/K_0N_0)^2 \le4N^2/M$ tiles. Thus, recalling the operation count of a single tile from Section 2.3.1, the whole product takes $O\left(\frac{N^2}{M}L \cdot10^L\right)$ operations, that is, $O(L \cdot10^L/M)$ per entry of $XY$. One can calculate that this is $N^2\operatorname{poly}(L)$ in total when $L = 10m$, but about $N^{2}D^{0.2}$ for our choice $L = 19m$. We nonetheless choose $L = 19m$, and in Section 2.4.3 below we explain why this is needed for our sparsity arguments.

Extracting only the desired entries

In this subsection, we present the main new idea of this paper. Our goal is now to compute only the wanted entries of $XY$, i.e., those at the positions in $W$, which is a set of size $|W| \le N^2/\sqrt{D}$. We will modify the algorithm of Section 2.3 in two ways. First, we will observe that the encoding (step (2)) of each input array depends only on that array and not the other, and that each input array is shared by many tiles (defined in Section 2.3.4), so we can compute each encoding only once and share it among those tiles (Section 2.4.1). Second, we will prune the recursion tree for the multiplication and decoding steps (steps (3) and (4)), by skipping a recursive call whenever none of its outputs leads to an output entry from $W$ (Section 2.4.2). Finally we will analyze the algorithm by counting the leaves that are reached (Section 2.4.3), and adding up the costs (Section 2.4.4).

Sharing the encoding

The two numbers $\Phi_{\tau}(a)$ and $\Psi_{\tau}(b)$ multiplied at a leaf $\tau$ depend on $a$ alone and on $b$ alone, respectively, so we will compute them once per input array rather than once per tile. The encoding of $a$ is the array of the $10^L$ numbers $\Phi_{\tau}(a)$, indexed by the leaves $\tau$. This is the array that the encoding steps of Full compute from $a$. To compute it, we run Full with $b$ left out and without step (4): by Lemma 7, the left number at the leaf $\tau$ is $\Phi_{\tau}(a)$, and the leaf stores it. This takes $O(10^L)$ operations, the cost of the encoding in Section 2.3.1. The encoding of $b$, the array of the numbers $\Psi_{\tau}(b)$, is computed in the same way.

We compute the encodings of the input arrays of all row bands and all column bands (defined in Section 2.3.4), $2N/K_0N_0$ arrays of $10^L$ numbers in all, and every tile reads its two encodings from these. It is important that these encodings are shared, since an encoding has $10^L=(9+1)^L$ numbers, more than the $M=\binom{L}{m}9^{L-m}$ entries of the $K$ products of a tile, so encoding every tile separately would cost more than one operation per entry of the entire matrix product (even before focusing on $W$). Since each encoding is computed only once per band, the $2N/K_0N_0$ encodings are used by all $(N/K_0N_0)^2$ tiles. Each one costs $O(10^L)$, independently of $W$. The assumption $N \ge D^{18}$ makes this cost negligible (as we will compute in Section 2.4.4).

Skipping the calls that are not needed

Once the encodings are computed as above, we can skip that phase of Full. Instead, each leaf reads its two numbers from the encodings, and each call only performs the decoding of step (4). We now modify Full so that each call computes only the outputs that are wanted from it, i.e., the outputs that lead to entries in $W$. We will implement this by passing to each call the set $S$ of the outputs that are wanted from it.

More precisely, the outputs of the vertex $\tau_1\cdots\tau_k$ (defined in Section 2.3.1) are indexed by output strings of length $L-k$, and $S$ is a set of such strings. If $U$ is the set passed to the root (in the final algorithm, this will be the set of output strings of the wanted entries of the tile), we will see in Lemma 10 that $S$ will consist of the suffixes $w_{k+1}\cdots w_L$ of the strings $w \in U$ to which the vertex contributes (as defined in Section 2.3.2). Sets of output strings have slices like arrays (see Section 2.3.1): for a set $S$ of output strings and an output variable $z$, the slice $S_z$ is the set of output strings $w'$ with $w'z' \in S$, and an array on $S$ assigns an integer to each string of $S$.

Our modified version of the recursive algorithm Full, which we call Pruned, is as follows. $\operatorname{Pruned}_{\tau_1\cdots\tau_k}(S)$, for a set $S$ of output strings of length $L-k$, returns the array returned by the vertex $\tau_1\cdots\tau_k$ of $\operatorname{Full}$, restricted to the strings in $S$:

1. If $k=L$, then look up $\Phi_\tau(a)$ and $\Psi_\tau(b)$ from the precomputed encodings (Section 2.4.1) and return $\Phi_\tau(a)\Psi_\tau(b)$ for the leaf $\tau=\tau_1\cdots\tau_L$.

2. For each term $\lambda$ of Schönhage’s identity, let $S_\lambda$ be the union of the slices $S_z$ over the $z$ to which $\lambda$ contributes, i.e., $S_{P_{ij}}:=S_{z_{ij}}\cup S_{z_0}$ and $S_{P_0}:=S_{z_0}$.

3. For each term $\lambda$ of Schönhage’s identity with $S_\lambda\ne\varnothing$, compute $C_\lambda:=\operatorname{Pruned}_{\tau_1\cdots\tau_k\lambda}(S_\lambda)$, an array on $S_\lambda$.

4. Return the array $c$ on $S$ whose slices are $c_{z_{ij}}:=C_{P_{ij}}$ and $c_{z_0}:=\sum_\lambda C_\lambda$, restricted to $S_{z_{ij}}$ and to $S_{z_0}$. (These restrictions make sense, since $S_{z_{ij}}\subseteq S_{P_{ij}}$ and $S_{z_0}\subseteq S_\lambda$ for every $\lambda$, by step (2).)

Let us confirm which leaves $\operatorname{Pruned}$ visits. Recall from (2) that $\operatorname{Mult}(a,b)[w]$ is the sum of $\Phi_\tau(a)\Psi_\tau(b)$ over the leaves contributing to $w$ (as defined in Section 2.3.2), i.e., those that choose $P_{ij}$ at every level where $w$ has $z_{ij}$, and may choose any term at levels where $w$ has $z_0$. For a set $U$ of output strings of length $L$, we write $\operatorname{Leaves}(U)$ for the set of leaves that contribute to some output string of $U$.

Lemma 10. Let $a$ and $b$ be input arrays, and let $U$ be a set of output strings of length $L$. Called at the root, $\operatorname{Pruned}(U)$ returns $\operatorname{Mult}(a,b)$ restricted to $U$. The leaves it visits are exactly those of $\operatorname{Leaves}(U)$, and the sets passed to the calls it makes have total size at most $(L+1)|\operatorname{Leaves}(U)|$.

Proof. The values are correct by induction on $L-k$: each $C_\lambda$ is the array of the corresponding child of $\operatorname{Full}$ restricted to $S_\lambda$, so step (4) is step (4) of $\operatorname{Full}$ restricted to $S$.

Recall that $S_\lambda$ consists of the last $L-k-1$ variables of each string in $S$ such that $\lambda$ contributes to the first variable of that string. Thus, an induction on $k$ shows that the set passed to the vertex $\tau_1\cdots\tau_k$ consists of the suffixes $w_{k+1}\cdots w_L$ of the strings $w$ of $U$ to which the vertex contributes, and furthermore that the vertex is called if and only if this set is nonempty. This means a leaf is visited if and only if it contributes to some string of $U$.

For the total size, map each suffix in the set passed to a vertex $\tau_1\cdots\tau_k$ to the leaf that extends $\tau_1\cdots\tau_k$ by picking $P_{ij}$ at the levels where the suffix has $z_{ij}$ and by picking $P_0$ where it has $z_0$. By the previous paragraph, the vertex contributes to some string $w$ of $U$ with this suffix, and then so does this leaf, so the leaf is in $\operatorname{Leaves}(U)$. The leaf determines the vertex and the suffix. This means the sets passed to the vertices at each of the $L+1$ depths have total size at most $|\operatorname{Leaves}(U)|$.

In other words, each call is made at most once, no matter how many strings of $U$ depend on it, and the running time depends on the size of the union of the sets of leaves contributing to the strings of $U$, rather than the sum of their sizes. It remains to bound this union.

Few leaves contribute to a sparse set of entries

In our final algorithm, which we present below in Section 2.4.4, the set $U$ passed to $\operatorname{Pruned}$ will consist only of the output strings we have been focusing on, i.e., output strings of length $L$ whose inner sets (defined in Section 2.3.3) have exactly $m$ elements, so that they index the entries of the products $X_QY_Q$ (Lemma 9). There are $M=\binom{L}{m}9^{L-m}=K N^2$ such output strings (with $M$ as defined in Section 2.3.4), and in this subsection all output strings are of this kind.

Our focus in this subsection is on counting how many leaves the $\operatorname{Pruned}$ algorithm visits. A single output string has $10^m=D^{\log_4 10}\approx D^{1.66}$ leaves contributing to it. This may appear problematic, since it is more than the $D$ multiplications that computing that output entry directly, as an inner product, would take. Because of this, it is important that we do not simply count the leaves of each output string of $U$ separately. Instead, we will find that most leaves contributing to an output string also contribute to many other output strings, and there are few such leaves in total.

Let us begin by classifying leaves by how many levels they choose $P_0$ at (as opposed to one of the $P_{ij}$s). As we mentioned above, the leaves contributing to an output string with inner set $Q$ are the $10^m$ leaves that choose any term at each of the $m$ levels of $Q$, but $P_{ij}$ at every level outside $Q$ where the output string has $z_{ij}$. We call the one that chooses $P_0$ at every level of $Q$ the private leaf of the output string. The private leaf contributes to no other output string, since an output string with a different inner set requires a term other than $P_0$ at some level of $Q$, and an output string with the same inner set but a different outer part (a different $z_{ij}$ at some level outside $Q$) requires a different term there.

To describe how close a leaf is to being a private leaf, we define the order of a leaf as $m$ minus the number of levels at which it chooses $P_0$. A leaf contributing to an output string chooses $P_0$ only at levels of the output string’s inner set, so its order is at least 0. Since $P_0$ is the only term that serves the inner product alone, this is the classification described in the technique overview: the order of a leaf contributing to an output string is the number of levels of the string’s inner set at which the leaf takes one of the terms $P_{ij}$ shared with the outer product. The only leaf of order 0 contributing to an output string is its private leaf. Beyond this, exactly

leaves of order $d$ contribute to an output string, namely its private leaf with $P_0$ replaced by one of the nine other terms at $d$ of the levels of $Q$ (Figure 6).

The total number of leaves of order $d$ in the entire recursion tree is

since we can choose the $m-d$ levels with $P_0$ and one of the nine other terms at every other level. In particular, $\beta_0 = \binom{L}{m}9^{L-m} = M$ since the leaves of order 0 are exactly the private leaves, and there is one per output string. A leaf of order $d$ contributes to $\binom{L-m+d}{d}$ output strings, one for each inner set of size $m$ that contains the $m-d$ levels at which it chooses $P_0$. Indeed, $M\alpha_d = \beta_d\binom{L-m+d}{d}$, as both sides count the pairs of an output string and a leaf of order $d$ contributing to it. Finally, if $L = 19m$, then for $1 \le d \le m$ we have

As in the technique overview, leaves of higher order are more numerous per output string, but each of them is shared by more output strings, and by (5) there are fewer of them in all. Below, we count all the leaves of order $d > m/9$, whether or not they contribute to a string of $U$, and this count will be small enough since the $\beta_d$ decrease geometrically. By the ratio in (5), this requires $m$ to be well below $L/10$: at $m \approx L/10$, the choice that computes the whole product cheaply (Section 2.3.4), the ratio would be close to 1 for small $d$, so the $\beta_d$ would stay about $M$ instead of decreasing, and the count would give no saving. A small $m$ also makes $D = 4^m$ small compared with $N_0 = 3^{18m}$ (about $N_0^{0.07}$), and since a tile must fit into the product, $N$ must be a large power of $D$. Together with the cost of the encodings (Section 2.4.4), this is why Theorem 5 assumes $N \ge D^{18}$.

From now on, we fix $L=19m$, so that $N_0=3^{18m}$ and $K=\binom{19m}{m}$. Recall that $D=4^m$, so $2^m=\sqrt{D}$ and $2^{-m/9}=D^{-1/18}$.

Lemma 11. For every set $U$ of output strings whose inner sets have exactly $m$ elements,

and if $L=19m$, then $\lvert\operatorname{Leaves}(U)\rvert\le2^{-m/9}(2^m\lvert U\rvert+2M)$.

Proof. We bound the number of leaves of order $d$ in $\operatorname{Leaves}(U)$ in two ways: by $\lvert U\rvert\alpha_d$, charging each output string of $U$ for its leaves of order $d$, and by $\beta_d$, the number of all leaves of order $d$. Since every leaf of $\operatorname{Leaves}(U)$ has an order $0\le d\le m$, summing the smaller of the two over $d$ gives the first inequality. (Summed on their own, $\lvert U\rvert\alpha_d$ gives $\lvert U\rvert\cdot10^m$, the cost of computing each output string of $U$ by itself, and $\beta_d$ gives less than $2M$ by (5), the cost of computing all $M$ output strings.) Neither $\lvert U\rvert\alpha_d$ nor $\beta_d$ depends on which output strings are in $U$.

For the second inequality, let $L=19m$. We use $\lvert U\rvert\alpha_d$ for the orders $d\le m/9$ and $\beta_d$ for the orders $d>m/9$ (Figure 7). The orders $d>m/9$ contribute at most $\sum_{d>m/9}\beta_d\le M\sum_{d>m/9}2^{-d}<2\cdot2^{-m/9}M$ by (5). The orders $d\le m/9$ contribute at most $\lvert U\rvert\sum_{d\le m/9}\binom{m}{d}9^d$, and since $9^d=72^d8^{-d}\le72^{m/9}8^{-d}$ for $d\le m/9$,

since $72\cdot(9/8)^9<208$. Adding the two parts gives the second inequality.

Since $2^m=\sqrt{D}$ and $2^{-m/9}=D^{-1/18}$, the second inequality of Lemma 11 says that

Consider a tile in which $|U|=M/\sqrt{D}$ entries are wanted. (This is the density of Theorem 5, which allows $N^2/\sqrt{D}$ wanted entries out of $N^2$.) Then $\sqrt{D}|U|=M$, and the bound says that only $O(M/D^{1/18})$ leaves contribute to the wanted entries of the tile. This is a factor $D^{1/18}$ below the number $M$ of output entries of the tile, and far below the $|U|\cdot D=M\sqrt{D}$ multiplications that computing each wanted entry directly, as an inner product, would take. Next in Section 2.4.4 we sum the number of leaves over all the tiles, and we will use the fact that the bound is linear in $|U|$, meaning the total is the same no matter how the wanted entries are spread over the tiles.

Proof of Theorem 5

We finally put everything together to state and analyze our algorithm. As before, set $D=4^m$ and $L=19m$.

Two consequences of $N \ge D^{18}$. We use the assumption $N \ge D^{18}$ twice, each time to check that $N$ is large enough compared with a quantity that grows exponentially in $m$. Since $D^{18}=4^{18m}$, we compare $m$-th powers by comparing their bases.

First, the tiling of Section 2.3.4 requires that $N \ge K_0N_0$ so that it does not lose too much when padding, since it pads $N$ to a multiple of $K_0N_0$. This holds: using the general bound $\binom{L}{m}\le(eL/m)^m$, we have $K=\binom{19m}{m}\le(19e)^m$, so $KN_0\le(19e\cdot3^{18})^m\le4^{18m}$, as $19e\cdot3^{18}<4^{18}$, and hence $K_0N_0\le KN_0\le D^{18}\le N$. This means the padding at most doubles $N$, and we continue to write $N$ for the padded size.

Second, the encoding phase of the algorithm (Section 2.4.1) has a cost that does not depend on $W$ at all, and we need to ensure it is not too large. Recall that we perform $2N/K_0N_0$ encodings, which use $O(10^L)$ operations each. To afford them within our budget of $O(2^{-m/9}N^2)$ operations, we need

This also holds: $10^{19}\cdot2^{1/9}/3^{18}<4^{18}$, so $10^{L}\cdot2^{m/9}/N_{0}\le D^{18}\le N$.

The algorithm. We compute the encodings of all row bands and column bands (Section 2.4.1). Each wanted position $(I,J)\in W$ lies in one tile $T$ and is indexed by one output string of $T$, the string whose inner set is the subset $Q$ of the block product containing $(I,J)$ and whose row and column are those of $(I,J)$ within that block product (Sections 2.3.3 and 2.3.4). For each tile $T$, let $W_T$ be the set of the output strings of the wanted positions in $T$. Then, for each tile $T$ with $W_T\ne\emptyset$, we run $\operatorname{Pruned}(W_T)$, and report, for each $(I,J)\in W$, the value at its output string. This is correct by Lemmas 9 and 10.

The cost. By (6), the encodings cost $(2N/K_0N_0)\cdot O(10^L)=O(2^{-m/9}N^2)$ operations. Each position of $W$ lands in one set $W_T$, so $\sum_T|W_T|\le|W|\le N^2/\sqrt{D}=2^{-m}N^2$, and there are at most $4N^2/M$ tiles, so $\sum_T M\le4N^2$. A call of $\operatorname{Pruned}$ with a set $S$ takes $O(L|S|)$ operations, since each of its steps handles each string of $S$ a constant number of times and a string has length $L$ (we store the sets as sorted lists, so that a slice is a segment and a union of slices is a merge). By Lemma 10, the total number of operations per leaf of $\operatorname{Leaves}(W_T)$ is thus $O(L^2)=O(\log^2D)$. Thus, by Lemma 11 summed over the tiles, the pruned recursions use

operations, no matter how $W$ is distributed among the tiles. The remaining minutiae, which list the $K_0^2$ subsets, locate the output strings, and sort the sets $W_T$, take $O(L)$ operations per subset, per tile, and per position of $W$, which is $O(L2^{-m/9}N^2)$ in all (as $K\le N$, there are at most $4N^2/N_0$ tiles, $|W|\le2^{-m}N^2$, and $N,N_0\ge2^{m/9}$). Finally, since $L=19\log_4D$ and $2^{-m/9}=D^{-1/18}$, the total operation count is $O(N^2\log^2D/D^{1/18})$, as desired.

Word size. Every number the algorithm handles is an integer of magnitude $N^{O(1)}$. The input entries are of this size, by assumption. Every number in an encoding, or formed on the way to one, is a $\pm1$ combination of at most $7^L$ input entries, as step (2) only adds and subtracts slices, and every value computed by $\operatorname{Pruned}$ is a sum of at most $10^L$ products of two encoded numbers (Lemma 7), where $7^L\le10^L\le N^2$ by (6). Indices and list positions are at most the number of words in use, also $N^{O(1)}$. Hence every integer has $O(\log N)$ bits, as desired. □

Remark 12 (comparison with related techniques). Prior work on matrix multiplication typically counts only multiplications, since they are usually a dominant cost for the algorithm, whereas in our algorithm, the additions need to be dealt with carefully. Our encodings (Section 2.4.1) are Yates’ algorithm for Kronecker powers [168, 93], as is standard in recursive matrix multiplication algorithms. The pruned recursion of Section 2.4.2 is an analog of FFT pruning [128, 142], of trimmed Möbius inversion [40], and, transposed, of sparse evaluation of Kronecker-power circuits [22], which builds on [167]. The new part along these lines in this paper is the count of Section 2.4.3, which shows that the pruned recursion visits few leaves.

Exact Triangle reduces to computing certain entries of a thin matrix product

This section derives the faster deterministic algorithms for Exact Triangle (given a tripartite graph with integer weights on its edges, is there a triangle whose three weights sum to zero?), 3SUM, and APSP from our new matrix theorem. Theorem 5 already suffices for this, and its stronger form in Section 4, Corollary 26, gives the running times stated in the introduction. Later, in Section 5, we will similarly derive our other algorithmic results.

Almost nothing in this section is new. We will apply known reductions from fine-grained complexity theory, due to Pătraşcu [135], Kopelowitz, Pettie, and Porat [116], Vassilevska Williams and Williams [158, 157, 159], Vassilevska Williams and Xu [160], Chan and He [51], Chan and Xu [71], and Fischer, Kaliciak, and Polak [84]. We go into some detail rather than just citing these known results for two main reasons. First, the theorem statements in the literature do not quite match the parameter regime of the matrix theorem, and obtaining exactly the matrix shapes that our algorithm needs requires slight modifications of the known constructions. Second, for one step, the reduction from Exact Triangle to the lopsided triangle problem, the only known reduction that achieves the parameters we need is randomized [160], and the known deterministic one [71] targets the balanced problem. We will make a small modification to achieve a deterministic reduction, using a known tool [84]. The reductions from 3SUM and APSP to Exact Triangle are deterministic already, and we cite them as they are. Thus, all the steps are deterministic, and so we ultimately design deterministic algorithms for all these problems. In Figure 8 we give an overview of the relevant reductions.

We remind the reader of some basics from fine-grained complexity that we will use here. A reduction from a problem $A$ to a problem $B$ is an algorithm for $A$ that builds instances of $B$, solves them using an oracle (an algorithm for $B$ that it treats as a constant time black box), and does some extra work. We will typically state three things about it: how many instances it builds, of what size, and how much extra time it uses. Composed with an algorithm for $B$, it is an algorithm for $A$, whose running time is the extra time plus the total time to solve $B$ on each instance created. Such a reduction can be interpreted as conferring hardness (if $A$ is hard then so is $B$) or as an algorithmic tool (a faster algorithm for $B$ gives one for $A$). Here we will use the second interpretation to design our algorithms.

The particular types of reductions we focus on are fine-grained reductions. Here, for two problems $A$ and $B$ and two time bounds $a(n)$, $b(n)$, an $(a(n),b(n))$-fine-grained reduction from $A$ to $B$ is an algorithm $R$ that solves any $n$-sized instance of $A$ by making oracle calls to instances of $B$ of some sizes $n_1,\ldots,n_t$ where for every $\delta> 0$ there is an $\varepsilon> 0$ s.t. $\sum_{i=1}^{t}(b(n_i))^{1-\delta} \leq(a(n))^{1-\varepsilon}$ and so that $R$ solves $A$ in $O(a(n)^{1-\varepsilon})$ time. In other words, if $B$ had an $O(b(N)^{1-\delta})$ time algorithm, then replacing the oracle calls by calls to the algorithm, due to the inequality above, we would get an $O(a(n)^{1-\varepsilon})$ time algorithm for $A$.

Section 3.1 restates Theorem 5 and Corollary 26 in the language of graphs. Section 3.2 reproves the reduction from Exact Triangle to the lopsided problem in a deterministic way. It hashes the weights modulo a small prime (in place of the random hash functions of earlier reductions: almost-linear hashing in the reductions from 3SUM [135, 116], and a random map into a large prime field in [160]), chooses the prime deterministically by counting false positives, as in the 3SUM reduction of Fischer, Kaliciak, and Polak [84], and builds the instances and searches for witnesses as in Chan and Xu [71]. The false positive counting in [84] is accomplished using the FFT; here, similar to prior work (e.g., [32]), for Exact Triangle we compute the count by multiplying two matrices whose entries are polynomials of low degree which can be done efficiently using fast matrix multiplication. Section 3.3 combines this reduction with Theorem 5 and Corollary 26 to solve Exact Triangle in truly subcubic time. Section 3.4 restates the reductions from 3SUM and APSP to Exact Triangle [158, 157, 159, 51] with their overheads, and applies them.

The Lopsided All-Edges Sparse Triangle problem

We define a lopsided (rectangular) version of the well-studied problem All-Edges Sparse Triangle. When viewed in graph terms, Theorem 5 counts the triangles through prescribed pairs of vertices in a tripartite graph whose middle part is small, and thus in particular can be used to solve All-Edges Sparse Triangle in such graphs.

Definition 13 (Lopsided All-Edges Sparse Triangle, Lop-AE-SparseTri($n$, $D$)). We are given an unweighted undirected tripartite graph with two parts $A$ and $B$ of $n$ vertices each and a middle part $M$ of at most $D$ vertices, with arbitrary edges in $M \times A$ and $M \times B$. Let $W \subseteq A \times B$ be the set of edges between $A$ and $B$. Decide, for every edge $(a,b) \in W$, whether it lies in a triangle with some vertex of $M$, that is, whether $a$ and $b$ have a common neighbor in $M$.

In this notation, the usual All-Edges Sparse Triangle problem is the balanced case $D = n$ (after the standard reduction to tripartite graphs), except that one asks about every edge of the graph, not only those between $A$ and $B$, and that the running time is measured in terms of the number $m$ of edges. We think of $W$ as being polynomially smaller than $n^2$ and $D$ as polynomially smaller than $n$, so the instance is sparse in this sense. The brute-force algorithm solves Lop-AE-SparseTri($n$, $D$) in $O(|W|D)$ time. If $D \le n^{0.321}$, then, since $0.321 \le\alpha$ [161],⁷ fast rectangular matrix multiplication can solve the problem in $n^{2+o(1)}$ time even when $|W| = n^2$. We will beat both $|W|D$ and $n^2$ for $D = n^\varepsilon$ and $|W| = n^2/\sqrt{D}$ with a constant $\varepsilon> 0$, and this will suffice to refute both the 3SUM and the APSP hypotheses.

Some of the known reductions, such as those from the real-valued problems in Section 5.2, reduce to the counting version of All-Edges Sparse Triangle, which asks for the number of triangles through each edge; we define the lopsided version.

Definition 14 (Lopsided All-Edges Sparse Triangle Counting, #Lop-AE-SparseTri($n$, $D$)). This is the same problem as Lop-AE-SparseTri($n$, $D$), except that for every edge ($a$, $b$) $\in W$ we ask for the number of triangles it lies in, that is, the number of common neighbors of $a$ and $b$ in $M$.

In matrix language, let $X \in\{0,1\}^{n \times D}$ and $Y \in\{0,1\}^{D \times n}$ be the biadjacency matrices of the edges between $A$ and $M$ and between $M$ and $B$, that is, $X[a,v] = 1$ if and only if $a \in A$ and $v \in M$ are adjacent, and $Y[v,b] = 1$ if and only if $v \in M$ and $b \in B$ are adjacent (padded with zeros if $M$ has fewer than $D$ vertices). Then Lop-AE-SparseTri($n$, $D$) asks for the entries of the Boolean product of $X$ and $Y$ at the positions in $W$, and #Lop-AE-SparseTri($n$, $D$) asks for the corresponding entries of the product of $X$ and $Y over the integers. In the language of Section 2, $W$ is the set of wanted positions, and we call its elements the query pairs.

The counting version #Lop-AE-SparseTri($n$, $D$) is also equivalent, up to polylogarithmic factors, to the problem of Theorem 5 and Corollary 26, that is, to computing the wanted entries $(XY)[I,J]$, $(I,J) \in W$, of a thin matrix product, for matrices $X \in\mathbb{Z}^{n \times D}$ and $Y \in\mathbb{Z}^{D \times n}$ whose entries have absolute value $n^{O(1)}$.⁸

Besides asking for the wanted entries of a Boolean thin matrix product, the detection version Lop-AE-SparseTri($n$, $D$) can also be viewed as offline Set Disjointness (see, e.g., [116, 160]): the sets are the neighborhoods in $M$ of the vertices of $A$ and $B$, over a universe of at most $D$ elements, and a query pair ($a$, $b$) $\in W$ asks whether the sets of $a$ and $b$ are disjoint. The reductions below (mostly from prior work) produce at most $n^2/\sqrt{D}$ query pairs and sets of at most $\sqrt{D}$ elements, which is the regime of [116] and of [160], Corollary 3.12.

With the above discussion in mind, Theorem 5 gives the following.

Corollary 15. Let $D \ge4$ be a power of four with $n \ge D^{18}$, and consider an instance of #Lop-AE-SparseTri($n$, $D$) or of Lop-AE-SparseTri($n$, $D$) with $|W|$ query pairs. If $|W| \le n^2/\sqrt{D}$, then the instance can be solved deterministically in $O(n^2 \log^2 D/D^{1/18})$ time. In general, splitting $W$ into sets of at most $n^2/\sqrt{D}$ query pairs solves it deterministically in time

Proof. Apply Theorem 5 with $N = n$ to the two biadjacency matrices: the entries $(XY)[a,b]$, $(a,b) \in W$, are the numbers of common neighbors, and they are nonzero exactly for the query pairs that lie in a triangle. For larger $W$, apply the theorem to each of at most $\lceil|W|\sqrt{D}/n^2\rceil+ 1$ pieces.

The data structure of Section 4 improves the saving $D^{1/18}$ to $D^{0.063}$ (Corollary 26), and gives the following.

Corollary 16. Let $n \ge D^{18}$, and consider an instance of #Lop-AE-SparseTri($n$, $D$) or of Lop-AE-SparseTri($n$, $D$) with $|W|$ query pairs. It can be solved deterministically in $O(|W|D^{0.437} + n^2/D^{0.063})$ time.

Proof. This is Corollary 26 with $N = n$, applied to the two biadjacency matrices as above.

For $|W| = O(n^2/\sqrt{D})$, which is what the reductions below produce, the running time in Corollary 16 is $O(n^2/D^{0.063})$.

The reduction of Section 3.2 is stated for an arbitrary $D$, with the number of oracle calls and the additional time as functions of $D$ and a free parameter $g$, and Section 3.3 chooses $D$ and $g$.

A deterministic reduction from Exact Triangle to Lop-AE-SparseTri

Exact Triangle. An instance of Exact Triangle (also called Zero-Weight Triangle, and sometimes “Love Triangle” [94]) consists of a complete tripartite graph on vertex parts $A, B, C$ of $n$ vertices each, and an integer weight $w(e)$ on every edge, with $|w(e)| \le n^\nu$ for some constant $\nu\ge1$. Let $S(a,b,c) := w(a,b) + w(b,c) + w(a,c)$. A zero triangle is a triangle $(a,b,c) \in A \times B \times C$ with $S(a,b,c) = 0$, and the task is to decide whether one exists. As is standard, an instance on an arbitrary $n$-vertex graph reduces to this form by taking three copies of the vertex set and giving every missing edge the weight $3n^\nu+ 1$, which lies in no zero triangle (this replaces $\nu$ by $\nu+ 1$), and the search version, which asks the algorithm to find a zero triangle if one exists, is equivalent to the decision version via a simple self-reduction (the reduction below finds one anyway).

Theorem 17 (Exact Triangle to Lop-AE-SparseTri, deterministically). Let $16 \le D \le n$, and let $1 \le g \le\sqrt{D}$ be an integer. Exact Triangle on $n$ vertices per part with weights of absolute value at most $n^\nu$ reduces deterministically to at most $4ng$ instances of Lop-AE-SparseTri$(n,D)$, each with at most $n^2/\sqrt{D}$ query pairs, plus $O(\nu n^3 \log n/g + n^{\omega+o(1)}D^{3/2} + n^2Dg)$ additional time, where $\omega< 2.372$ is the matrix multiplication exponent [10]. The reduction is non-adaptive: it produces all the instances before it makes any oracle call, so they can be solved in any order, or all together.

The parameter $g$ trades the number of instances against the cost of the search for witnesses in the last step of the proof (the term $\nu n^3 \log n/g$); Section 3.3 balances the two. The oracle returns at most $n^2/\sqrt{D}$ answers per instance, and the time to read them is counted as part of the oracle calls, rather than as additional time.

Proof. The proof has three steps: we hash the weights modulo a prime, build the instances, and search for witnesses.

Hashing modulo a prime. We reduce the weights modulo a prime $p \in[\sqrt{D}/2,\sqrt{D})$, chosen deterministically. This size of $p$ fits the instances below: a piece of $C$ of up to $\sqrt{D}$ vertices, paired with the $p$ labels, gives a middle part of at most $D$ vertices, and the $n^2$ pairs $(a,b)$ fall into $p$ residue classes of about $n^2/\sqrt{D}$ pairs on average, which is the number of query pairs that an instance may have. A false positive of $p$ is a triple $(a,b,c)$ with $S(a,b,c) \ne0$ and $p \mid S(a,b,c)$; let $F(p)$ denote the number of false positives of $p$. For every prime $p$ in the range we count the triples with $S(a,b,c) = 0 \pmod p$. This number is the number of false positives $F(p)$ plus $Z_0$, the number of zero triangles, which does not depend on $p$. Each $F(p) + Z_0$ is computed by a standard trick for small weights [14, 172]: let $P[a,c] := x^{w(a,c)} \bmod p$ and $Q[c,b] := x^{w(b,c)} \bmod p$ be matrices over the ring $\mathbb{Z}[x]/(x^p - 1)$, so that the coefficient of $x^r$ in $(PQ)[a,b]$ is the number of $c \in C$ with $w(a,c) + w(b,c) \equiv r \pmod p$; then $F(p) + Z_0$ is the sum over the pairs $(a,b) \in A \times B$ of the coefficient of $x^{-w(a,b)} \bmod p$ in $(PQ)[a,b]$. Computing $PQ$ takes $n^{\omega+o(1)}$ ring operations, or $O(n^{\log_2 7})$ with Strassen’s algorithm [147], each of them $O(p^2)$ word operations, so $n^{\omega+o(1)}D^{3/2}$ time over the fewer than $\sqrt{D}$ primes in the range. (Any matrix multiplication algorithm fast enough to make this term negligible suffices. For our choice of $D$, Strassen’s algorithm already does, and keeps the reduction potentially practical.⁹

We select the prime with the smallest count, which is also the prime with the fewest false positives, and call it $p$. To bound $F(p)$, note that a triple with $S(a,b,c) \ne0$ is a false positive exactly of the primes in the range that divide $S(a,b,c)$. Since $0 < \lvert S(a,b,c)\rvert\le3n^\nu$ and these primes are at least $\sqrt{D}/2$, there are at most $\log_{\sqrt{D}/2}(3n^\nu)$ of them. Hence the numbers of false positives of all the primes in the range add up to at most $n^3\log_{\sqrt{D}/2}(3n^\nu)$. By the prime number theorem there are $\Omega(\sqrt{D}/\log D)$ primes in the range, and $F(p)$ is at most the average over them, so, as in [84], $F(p)=O(n^3\log(3n^\nu)/\sqrt{D})=O(\nu n^3\log n/\sqrt{D})$, since $\log D=O(\log(\sqrt{D}/2))$ for $D\ge16$.

The instances. As in [71], let $s:=\lfloor\sqrt{D}\rfloor$ and split $C$ into pieces $C_1,\ldots,C_h$ of at most $\lceil s/g\rceil$ vertices each, so that $h\le\lceil ng/s\rceil$. For $\varrho\in\mathbb{Z}_p$ let $W_\varrho$ be the set of edges $(a,b)\in A\times B$ with $w(a,b)\equiv\varrho\pmod p$, and cut it into chunks of at most $n^2/\sqrt{D}$ query pairs. There are at most $p+\sqrt{D}\le2\sqrt{D}$ chunks in all. For each chunk $Q\subseteq W_\varrho$ and each piece $C_k$ form the instance of Lop-AE-SparseTri$(n,D)$ with query pairs $W:=Q$ and middle part $C_k\times\mathbb{Z}_p$, of size at most $sp\le D$, whose middle vertices are pairs (vertex, label) as in [71], with

Within the chunk, the condition $S(a,b,c)\equiv0\pmod p$ has become the equality $w(a,c)+\varrho\equiv-w(b,c)$ of a label of $(a,c)$ and a label of $(b,c)$, so a query pair $(a,b)\in Q$ has a common neighbor if and only if some $c\in C_k$ has $S(a,b,c)\equiv0\pmod p$. There are at most $2\sqrt{D}h\le2\sqrt{D}(ng/s+1)\le4ng$ instances, since $s\ge\sqrt{D}-1\ge3\sqrt{D}/4$ because $D\ge16$, and $\sqrt{D}\le\sqrt{n}\le2n/3$ because $D\le n$ and we may assume $n\ge3$, and they are produced non-adaptively, since $p$ is chosen before any of them. Writing them down costs $O(n^2Dg)$: the two bipartite graphs of an instance have $O(nD)$ entries, and the chunks are computed once and shared by the pieces.

Witnesses. For every query pair that the oracle accepts, scan the piece $C_k$ of its instance for a $c$ with $S(a,b,c)=0$, as in the exhaustive search of [71]. A zero triangle is always found, since its $c$ lies in some piece and makes the oracle accept its pair. We stop as soon as a zero triangle is found. A failed scan, of a piece $C_k$ for a pair $(a,b)$, contains a $c\in C_k$ with $S(a,b,c)\equiv0\pmod p$ but $S(a,b,c)\ne0$, that is, a false positive $(a,b,c)$ of $p$, and distinct scans contain distinct false positives, so there are at most $F(p)$ failed scans, and the scans cost $O((F(p)+1)\sqrt{D}/g)=O(\nu n^3\log n/g)$.

Remark 18 (Comparison with [71] and [160]). We emulate the instances and the witness search of Chan and Xu [71], Section 3. The difference is in how the pairs are grouped so that the condition on $S(a,b,c)$ becomes an equality of labels. For exact integer weights, $w(a,b)$ takes too many values to group by, so Chan and Xu group the pairs by a reference witness $k$ with known $S(a,b,k)$, which they find by a recursion on the bits of the weights, and use Fredman’s trick. Their labels are exact and there are no false positives. After hashing, $w(a,b)\bmod p$ takes only $p$ values, so we can group by it directly, in a single round of instances instead of a recursion. The labels take only $p$ values, so the middle part is small, but in exchange, this introduces false positives. Vassilevska Williams and Xu [160] also hash, but by a random map $w(u,v) \mapsto xw(u,v)+y_u-y_v$ into a large prime field, with $x$ and the vertex weights $y_v$ random, followed by a split of the field into intervals. With residues in place of their intervals, the block $C_k \times\{-j\}$ of our instance for $\varrho$ and $C_k$ is their graph $G_{-\varrho-j,j,\varrho}$, and our instance puts the $p$ graphs with the same $\varrho$ side by side. They use randomness to bound the degrees and the false positives through each edge, which their listing step needs [160], whereas the scans above need only the total number of false positives, controlled by the choice of $p$.

Exact Triangle in truly subcubic time

Theorem 19 (Exact Triangle). Let $\varepsilon' := 0.00175$ and $\varepsilon_T := 0.0017$. For every constant $\nu\ge1$, Exact Triangle on $n$ vertices per part with integer weights of absolute value at most $n^\nu$ can be solved by a deterministic algorithm in $O(n^{3-1/648}\log^2 n)$ time using Theorem 5 (via Corollary 15), and in $O(n^{3-\varepsilon'}\log n) \le O(n^{3-\varepsilon_T})$ time using Corollary 26 (via Corollary 16).

Proof. Assume $n$ is larger than a constant depending on $\nu$ (smaller instances are solved by brute force). In both cases we apply Theorem 17, whose instances have at most $n^2/\sqrt{D}$ query pairs each, and solve every instance by a corollary of Section 3.1. The values of $g$ below balance the cost of the instances against that of the scans (Remark 20).

By Theorem 5. Let $D$ be the largest power of four with $D \le n^{1/18}$, so that $D \ge n^{1/18}/4$, and let $g := \lceil D^{1/36}\rceil$. Solve every instance by Corollary 15, which applies because $D^{18} \le n$. The $O(nD^{1/36})$ instances cost $O(n^2\log^2 D/D^{1/18})$ each, so $O(n^3D^{-1/36}\log^2 n)$ in all, the scans cost $O(n^3D^{-1/36}\log n)$, the choice of $p$ costs $O(n^{\log_2 7}D^{3/2}) = O(n^{2.9})$ with Strassen’s algorithm, and building the instances costs $O(n^2D^{1.03})$. Finally $D^{-1/36} \le4^{1/36}n^{-1/648}$, so the time is $O(n^{3-1/648}\log^2 n)$.

By Corollary 26. Let $D := \lfloor n^{1/18}\rfloor$ and $g := \lceil D^{0.0315}\rceil$, and solve every instance by Corollary 16, which applies because $D^{18} \le n$. The instances cost $O(nD^{0.0315}) \cdot O(n^2/D^{0.063}) = O(n^3D^{-0.0315})$, the scans $O(n^3D^{-0.0315}\log n)$, the choice of $p$ costs $O(n^{\log_2 7}D^{3/2}) = O(n^{2.9})$ with Strassen’s algorithm, and building the instances costs $O(n^2D^{1.04})$. Finally $D \ge n^{1/18}/2$ gives $D^{-0.0315} \le2n^{-0.0315/18} = 2n^{-0.00175}$, so the time is $O(n^{3-\varepsilon'}\log n) \le O(n^{3-\varepsilon_T})$.

Remark 20. With $g = D^\eta$ and a saving $D^\gamma$ per instance, the instances cost $n^3D^{\eta-\gamma}$ and the scans $n^3D^{-\eta}$, up to logarithmic factors. They balance at $\eta= \gamma/2$, that is, at $\eta= 1/36$ for the saving $D^{1/18}$ of Theorem 5 and at $\eta= 0.0315$ for the saving $D^{0.063}$ of Corollary 26. Any improvement of these savings improves the final exponent $\eta/18$ directly, and the factor $1/18$ in it comes from the assumption $N \ge D^{18}$ of these theorems. Any other algorithm for Lop-AE-SparseTri$(n,D)$ can be plugged into Theorem 17 in the same way. With the straightforward algorithm ($O(n^2/g)$ per instance, since every vertex of $A$ has $O(\sqrt{D}/g)$ neighbors in the middle part) the reduction recovers the brute-force bound $n^3$ up to a logarithmic factor.

3SUM and APSP reduce to Exact Triangle

The reductions here are known and deterministic. We recall them in the form we use, and then apply Theorem 19, with its bound $O(n^{3-\varepsilon'}\log n)$, $\varepsilon' = 0.00175$, when it uses Corollary 26.

Theorem 21 (Known reductions). For every constant $\nu$: (a) (3SUM) $3SUM$ on $n$ integers of absolute value at most $n^\nu$ reduces deterministically, in $n^{3/2+o(1)}$ time, to $n^{1/2+o(1)}$ instances of Exact Triangle on $n^{1/2+o(1)}$ vertices per part with weights of absolute value $n^{O(1)}$ [51, 158].

(b) ((min, +)-product and APSP) If a deterministic algorithm solves Exact Triangle on $s$ vertices per part with weights of absolute value at most $cU$, for a suitable constant $c$, in time $T(s)$ with $T(s)/s$ nondecreasing, then the (min, +)-product of two $n \times n$ integer matrices with entries of absolute value at most $U$ can be computed deterministically in $O(n^2T(n^{1/3})\log^2 U)$ time, and APSP on directed $n$-vertex graphs with integer weights of absolute value at most $n^\nu$ and no negative cycles in $O(n^2T(n^{1/3})\log^3 n)$ time [157, 159, 158].

Using the above known reductions we solve both 3SUM and APSP polynomially faster:

Theorem 22 (3SUM and APSP). Let $\nu$ be a constant. Using Theorem 5 (through Theorem 19), deterministic algorithms solve $3SUM$ on $n$ integers of absolute value at most $n^\nu$ in $n^{2-1/1296+o(1)} \le O(n^{1.99923})$ time, and the (min, +)-product of two $n \times n$ integer matrices with entries of absolute value at most $n^\nu$, as well as APSP on directed $n$-vertex graphs with integer weights of absolute value at most $n^\nu$ and no negative cycles, in $\widetilde{O}(n^{3-1/1944}) \le O(n^{2.99949})$ time. Using Corollary 26 instead, the times are $n^{2-\varepsilon'/2+o(1)} \le O(n^{1.9992})$ and $\widetilde{O}(n^{3-\varepsilon'/3}) \le O(n^{2.99942})$.

Proof. Plug Theorem 19 into Theorem 21. With the bound $O(n^{3-1/648}\log^2 n)$ of Theorem 19, the cost for 3SUM is

and for the (min, +)-product and APSP, with $T(s)=O(s^{3-1/648}\log^2 s)$, for which $T(s)/s$ is nondecreasing, it is

With the bound $O(n^{3-\varepsilon'}\log n)$ instead, the cost for 3SUM is

and for the (min, +)-product and APSP, with $T(s)=O(s^{3-\varepsilon'}\log s)$, it is

Remark 23 (3XOR). The same approach as for 3SUM gives a deterministic truly subquadratic algorithm for 3XOR (given three lists of $n$ vectors in $\mathbb{F}_2^{O(\log n)}$, find one vector from each list such that the three XOR to zero [108, 77]), using linear hash functions over $\mathbb{F}_2$ whose rows are chosen one at a time from an $\varepsilon$-biased set [13] by the method of conditional expectations.

The matrix theorem in general: a data structure

In this section, we extend our algorithm from Section 2 in two directions. First, we turn it into a data structure. After preprocessing $X$ and $Y$ in less time than it takes to write down $XY$, the data structure answers a query for any single entry of $XY$, not known in advance, in time polynomially smaller than the $D$ operations that it would take to compute that entry as an inner product. Second, we let the parameters of our construction vary. This will give larger time savings than $D^{1/18}$, and trade-offs between the time savings, how thin the product must be, the query time, and the density of the set $W$ of wanted entries. We assume that the reader is familiar with Section 2, whose notation we use throughout, although we will begin with a brief reminder of the relevant notions.

In Section 4.1, we state the main results, Theorems 24 and 25, with explicit parameter examples in Table 2, as well as Corollary 26, the instance of Theorem 24 that our reductions in Sections 3 and 5 use. In Section 4.2, we will introduce the key new idea, the partial sums that our data structure stores, which we call boxes. In Section 4.3, we will then give the data structure in terms of the parameters (listed in Table 1), and finally we will set the parameters in Section 4.4.

symbol	meaning	in Section 2	in this section
$m$, $D = 4^m$	the number of inner levels, and the inner dimension of the product $XY$
$L$, $c = L/m$	the number of levels of the recursion, and its ratio to $m$	$L = 19m$	$L = cm$, $c > 10$
$N_0 = 3^{L-m}$	the number of rows of a row block, and of columns of a column block
$K = \binom{L}{m}$	the number of products $X_QY_Q$ of a tile, one for each set $Q$ of $m$ levels
$M = KN_0^2$	the number of output entries of a tile
$\alpha_d = \binom{m}{d}9^d$	the number of leaves of order $d$ contributing to one output entry
$\beta_d = \binom{L}{m-d}9^{L-m+d}$	the number of leaves of order $d$ in a tile	$\beta_d \le2^{-d}M$	$\beta_d \le\rho^dM$
$\rho= 9m/(L-m+1)$	the decay rate of the $\beta_d$: $\beta_d/\beta_{d-1} \le\rho$, see (7)	$\rho< \frac{1}{2}$	$\rho< 1$
$t$	the switching order: a query reads the leaves of order below $t$, and the preprocessing puts the others into boxes; in Section 2, the order at which the count of Lemma 11 splits	$t = m/9$	$t = \theta m$, $0 < \theta< 0.9$
$\gamma$	the preprocessing, or the whole computation when $W$ is given in advance, takes $O(N^2 \log^2 D/D^\gamma)$ time	$\gamma= 1/18$	see Table 2
$\varepsilon$	the results hold for $D \le N^\varepsilon$	$\varepsilon= 1/18$	$\varepsilon< R_c(\gamma) < 0.1204\ldots$
$q$	a query takes $O(D^q \log D)$ time	no queries	see Table 2
$\kappa$	the set $W$ of wanted entries has $|W| \le N^2/D^\kappa$	$\kappa= \frac{1}{2}$	see Table 2

ratio $c$	$D \le N^{\varepsilon}$ $\varepsilon$	\multicolumn{6}{c}{Theorem 24: the data structure $\gamma$ for query time $D^q$, $q =$}	\multicolumn{5}{c}{Theorem 25: the wanted entries $\gamma$ for $|W| \le N^2/D^\kappa$, $\kappa=$}
		0.10	0.25	0.43	0.50	0.75	0.90	0.10	0.25	0.50	0.75	1.00
40	0.029	0.0205	0.0608	0.1175	0.1429	0.2417	0.3095	0.0166	0.0473	0.1060	0.1720	0.2442
21	0.056	0.0112	0.0331	0.0640	0.0778	0.1316	0.1685	0.0099	0.0286	0.0653	0.1074	0.1546
19	0.062	0.0097	0.0287	0.0555	0.0675	0.1142	0.1463	0.0087	0.0253	0.0578	0.0954	0.1379
15	0.079	0.0061	0.0183	0.0354	0.0430	0.0728	0.0932	0.0057	0.0168	0.0389	0.0646	0.0941
12	0.099	0.0028	0.0083	0.0160	0.0195	0.0330	0.0423	0.0027	0.0080	0.0186	0.0312	0.0459
10.5	0.114	0.0007	0.0022	0.0043	0.0052	0.0089	0.0114	0.0007	0.0022	0.0052	0.0087	0.0129

Notions from Section 2. We keep the recursion, the tiling, and the encodings of Sections 2.3 and 2.4.1, and we use the terminology of Sections 2.3.2 and 2.4.3 as defined there: the leaves contributing to an output string, its private leaf, the order of a leaf, and the counts $\alpha_d$ and $\beta_d$. For convenience, we briefly recall these notions here.

• A leaf $\tau= \tau_1 \cdots\tau_L$ is a string of $L$ terms of Schönhage’s identity, where $\tau_\ell$ is the term chosen at level $\ell$ of the recursion (Section 2.3.1). Throughout this section, we compare levels by their indices: we say that level $\ell$ is lower than level $\ell'$ if $\ell< \ell'$, so that the lowest level is the first level of the recursion.

• An output string $w$ with inner set $Q$ has the variable $z_0$ at each of the $m$ levels of $Q$, and one of the variables $z_{ij}$ at every other level. It indexes the entry $(X_QY_Q)[w]$ of the product $X_QY_Q$ (Section 2.3.3).

• A leaf contributes to $w$ if it chooses $P_{ij}$ at every level where $w$ has $z_{ij}$, while it may choose any of the ten terms at the levels of $Q$, where $w$ has $z_0$ (Section 2.3.2).

• The private leaf of $w$ is the leaf contributing to $w$ that chooses $P_0$ at every level of $Q$ (Section 2.4.3).

• The order of a leaf is $m$ minus the number of levels at which it chooses $P_0$. For a leaf contributing to $w$, this is the number of levels of $Q$ at which it chooses a term other than $P_0$ (Section 2.4.3).

• A tile consists of a band of $K_0 = \lfloor\sqrt{K}\rfloor$ row blocks of $X$ together with a band of $K_0$ column blocks of $Y$ (Section 2.3.4), and the encoding of a band is the array of the $10^L$ numbers $\Phi_\tau(a)$, or $\Psi_\tau(b)$, computed from its input array (Section 2.4.1). We will recall both in more detail in Section 4.3.

As in Section 2.4.3, all output strings in this section have inner sets of exactly $m$ elements, so that they index the $M = KN_0^2$ entries of the $K$ products $X_QY_Q$ of a tile, and we will use that $M = \beta_0 \le10^L$. For a leaf $\tau$ of a tile with input arrays $a$ and $b$, we call $\Phi_\tau(a)\Psi_\tau(b)$ the *product at $\tau$. We read its two factors from the encodings of the tile (Section 2.4.1), and by (2) and Lemma 9, $(X_QY_Q)[w]$ is the sum of the products at the leaves contributing to $w$. The one change from Section 2 is that $L$ is now any integer with $L \ge10m$, and we call $c = L/m$ the ratio.

Technique overview. The pruned recursion of Section 2.4 needs to know $W$ in advance, since it visits exactly the leaves that the wanted entries need. When the queries arrive one at a time, we must instead prepare for every possible output entry.

Recall that our preprocessing must take less time than it takes to write down $XY$. This means it stores fewer than one number per entry of $XY$, so each number that it stores must be useful for many different entries. Meanwhile, a query must take far fewer than the $D$ operations that it would take to compute the entry directly as an inner product. In Section 2, the running time had three parts: the leaves of small order, the leaves of large order, and the encodings (Lemma 11 and Section 2.4.4). Our data structure will divide these parts between the queries and the preprocessing. Each output entry uses only a few of the leaves of small order, so a query may just read them one by one. By contrast, each output entry may use many leaves of large order, but there are not too many such leaves in total, and each of them is shared by many output entries. Thus, in preprocessing we will add them up in advance, in groups which can be quickly summed during queries. The preprocessing also computes the encodings, once.

The $10^m$ leaves contributing to an output string form a cube: at each level of the inner set $Q$ of the output string, they choose any of the ten terms, and at each other level, they agree with the private leaf of the output string. Every one of these leaves, other than the private leaf, also contributes to other output entries. We will therefore add up parts of every cube in advance, into what we call boxes (which are themselves smaller cubes), so that a query adds up a few box values instead of $10^m$ products. For this we need a rule that cuts every cube into few boxes, with each of its leaves of large order in exactly one of them; we will prove this for our rule in Lemma 27. We will make sure that the definition of a box does not refer to the output entry, so that the same box serves many output entries, which keeps the total number of boxes small, as we will show in Lemma 29.

We also need to compute all the boxes of a tile quickly in preprocessing, and here we will use dynamic programming. Since a box is a cube, it has a fixed term at some levels, and its other levels are free, meaning all ten terms are allowed there. Fixing each of the ten terms in turn at the highest free level splits the box into ten smaller boxes, and its value is the sum of their values. Thus, in Lemma 29, we will compute the boxes from the smallest up, each from the values of smaller boxes computed before.

The order $t$ at which a query switches from leaves to boxes, which we call the switching order, trades off the query time (about $\alpha^t$) against the preprocessing time (about $\rho^t$ boxes per entry of the product, where $\rho< 1$ is the rate at which the $\beta_d$ decay). The ratio $c = L/m$ (which we fixed to be 19 in Section 2, but we now allow to vary) trades off the time savings against how thin the product must be. While a larger $c$ makes $\rho$ smaller, it also makes the encodings larger relative to $D$, so that computing them stays within the time bound only when $N$ is a larger power of $D$. With $c = 21$ and $t = m/9$, we will get the instance $N \ge D^{18}$ of Corollary 26, and pushing $c$ down toward 10 will give the upper limit for $D$ of our technique, namely $D \le N^{0.1204}$ (Corollary 31).

The results

In this section, we design a data structure with the following guarantee. Let $\varepsilon^* := \ln4/(5 \ln10) = 0.1204\ldots$.

Theorem 24. For every $\varepsilon< \varepsilon^*$, and every $q > 0$, there is a $\gamma> 0$ such that the following holds. Given as input matrices $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$, where $2 \le D \le N^\varepsilon$, whose entries are integers of absolute value at most $N^{O(1)}$, we can preprocess them deterministically in $O(N^2 \log^2 D/D^\gamma)$ time and space, after which any single entry $(XY)[I,J]$ can be computed deterministically in $O(D^q \log D)$ time.

The constant $\varepsilon^*$ gives the upper limit for $D$ of our technique, and we will explain where it comes from in Section 4.4 below (after the proof of Corollary 31). To achieve the general form of the matrix theorem, where we know the positions in advance as in Section 2, we simply ask the data structure one query per position. Setting the parameters to balance the preprocessing time and the total time for all the queries gives the following.

Theorem 25. For every $\varepsilon< \varepsilon^*$ and every $\kappa> 0$ there is a $\gamma> 0$ such that the following holds. Given as input matrices $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$, where $2 \le D \le N^\varepsilon$, whose entries are integers of absolute value at most $N^{O(1)}$, as well as a set $W$ of at most $N^2/D^\kappa$ positions of an $N \times N$ matrix, the entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in time $O(N^2 \log^2 D/D^\gamma)$.

The logarithmic factors in both theorems can be removed by halving $\gamma$ and applying Theorem 24 with $q/2$ in place of $q$; for $D = 1$ the bounds are trivial.

In Table 2, we list several parameter settings for Theorems 24 and 25. Theorem 5 is the entry $c = 19$, $q = 0.43$, applied to a set $W$ of $N^2/\sqrt{D}$ wanted entries (that is, $\kappa= \frac{1}{2}$). We next state explicitly the parameter setting that we use in the introduction and in our reductions in Sections 3 and 5. It is the entry $c = 21$, $q = 0.43$ of Table 2, whose switching order is $t = m/9$, as in Section 2, with the logarithmic factors absorbed into the exponents.

Corollary 26. Let $N \ge D^{18}$, and let $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$ have entries of absolute value at most $N^{O(1)}$. We can preprocess them deterministically in $O(N^2/D^{0.063})$ time and space, after which any single entry $(XY)[I,J]$ can be computed deterministically in $O(D^{0.437})$ time. Hence, for every set $W$ of positions of an $N \times N$ matrix, the entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in $O(|W|D^{0.437} + N^2/D^{0.063})$ time, which is $O(N^2/D^{0.063})$ whenever $|W| \le N^2/\sqrt{D}$.

Boxes

In this subsection, we define the boxes and prove two important facts about them. Intuitively, a box is a group of leaves that is shared by many output strings. In preprocessing we will add up the products at the leaves of each box, and then that value can help to answer many queries. In terms of the switching order $t$, we will ensure that the boxes of an output string partition its leaves of order at least $t$, and that there are only $\alpha_t$ of them. We will thus see, in Lemmas 27 and 28, that every output entry is the sum of a few products and a few values of boxes.

Cubes. A cube is a string $\pi= \pi_1 \cdots\pi_L$ of $L$ symbols, each of which is one of the ten terms of Schönhage’s identity or a star $*$. Its leaves are the leaves obtained by replacing each star by any of the ten terms, so a cube with $e$ stars has $10^e$ leaves, and a cube without stars is a single leaf.

The value of a cube $\pi$ is the sum of the products at its leaves,

Every output entry is the value of a cube. Indeed, recall that the entry $(X_QY_Q)[w]$ of an output string $w$ with inner set $Q$ is the sum of the products at the $10^m$ leaves contributing to $w$. These are the leaves of the cube of $w$, the cube $\pi^w$ that has a star at each level where $w$ has $z_0$ (the levels of $Q$) and the term $P_{ij}$ at each level where $w$ has $z_{ij}$, so that $(X_QY_Q)[w] = \operatorname{val}(\pi^w)$. In Figure 6, this cube is all the rows taken together. Figure 9 shows the cube of an output string, and the boxes defined next, in the recursion tree of a small example.

We will cut the cube of every output string into parts that are again cubes, in a way where each leaf of the cube lies in exactly one part, so that

and the value of a part will be shared by all the output strings that have it as a part. In Lemma 28 below, we will prove this equation for the parts that we will choose. There will be two kinds of parts. Each leaf of small order will be a part on its own (recall that a single leaf is a cube without stars), whereas the leaves of large order will be grouped together into larger parts, which we call boxes.

The definition of a box. The order of a leaf tells us how widely it is shared between output entries. Recall that the order is $m$ minus the number of levels at which the leaf chooses $P_0$, so that $\alpha_d = \binom{m}{d}9^d$ leaves of order $d$ contribute to an output string, and that a leaf of larger order contributes to more output strings. Throughout, we fix the switching order $t$, an integer with $0 \le t \le m$. The boxes will collect the leaves of order at least $t$, i.e., the leaves that choose $P_0$ at most $m-t$ times, and a query will read the leaves of order below $t$ one by one (Section 4.3).

Before giving the formal definition, we describe the idea informally and give some examples. Consider an output string $w$ with inner set $Q$, and a leaf $\tau$ of order at least $t$ contributing to $w$. To find the box of $\tau$, we look at the levels of $Q$ one at a time, from the highest down, and we stop as soon as we have seen $t$ levels at which $\tau$ chooses a term other than $P_0$. The box of $\tau$ then agrees with $\tau$ everywhere, except that it has a star at every level of $Q$ below the level where we stopped. For instance, let $m=4$ and $t=2$ as in Figure 10, and suppose that $\tau$ chooses the terms $P_{12}, P_0, P_{31}, P_{22}$ at the four levels $\ell_1 < \ell_2 < \ell_3 < \ell_4$ of $Q$. We see $P_{22}$ at level $\ell_4$ and $P_{31}$ at level $\ell_3$, and we stop there, so the box of $\tau$ is $* \ * \ P_{31}\ P_{22}$ (in the first row of the figure). This box contains $\tau$ along with the 99 other leaves that differ from $\tau$ only at the levels $\ell_1$ and $\ell_2$. Similarly, for a leaf that chooses the terms $P_0, P_{23}, P_0, P_{11}$, we stop at level $\ell_2$, so its box is $* \ P_{23}\ P_0\ P_{11}$ (in the second row of the figure). Since a box has stars at the levels that we did not reach, it contains many leaves. Moreover, although we described the box of $\tau$ in terms of $w$, the definition of a box below does not refer to an output string, and the same box serves many output strings.

A box is a cube $\pi=\pi_1\cdots\pi_L$ such that (i) at most $m-t$ of its symbols are $P_0$ or stars, and (ii) every star is at a lower level than every $P_0$. That is, the symbols $P_0$ and stars of a box, read in increasing order of the levels, are

and its other $L-f$ symbols are among the nine terms $P_{ij}$. A leaf of $\pi$ chooses $P_0$ at the $f-e$ levels where $\pi$ has $P_0$ and possibly at some of its $e$ star levels, and hence at most $f\le m-t @$Winvalid? times by (i), so a box only contains leaves of order at least $t$. By (ii), a box is obtained from a leaf of order at least $t$ by replacing its $e$ lowest symbols $P_0$ (rather than an arbitrary subset of them) by stars. Therefore, each leaf yields at most $m+1$ boxes, one for each value of $e$, which keeps their total number small, as we will show in Lemma 29. Boxes of this shape suffice: we will prove in Lemmas 27 and 28 that the leaves of order at least $t$ of an output string can be split among $\alpha_t$ of them. (We will see below that the boxes of an output string all have exactly $m-t$ symbols that are $P_0$ or stars. We nonetheless allow fewer in the definition of a box, since boxes with fewer such symbols will arise as intermediate values in the dynamic program of Lemma 29.)

It is important here that a star stands for all ten terms, $P_0$ included, so that one box collects leaves of several orders (from $m-f$ to $m-f+e$, for a box with $e$ stars and $f-e$ symbols $P_0$). This is needed for fast queries. If a star stood only for the nine terms $P_{ij}$, then every box would fix the set of levels at which its leaves choose $P_0$, so an output string would need a separate box for each subset of $Q$ of size at most $m-t$. That would be about $2^m=\sqrt{D}$ boxes when $t$ is small, so a query would take about $\sqrt{D}$ time, which gives no saving below $N^2$ when there are $N^2/\sqrt{D}$ queries.

The boxes of an output string. We now make precise the informal rule above for finding the box of a leaf, and describe all the boxes of an output string $w$ with inner set $Q$. Let $\tau$ be a leaf of order at least $t$ contributing to $w$, whose box we found above by going through the levels of $Q$ from the highest down (Figure 10). In symbols, let

(the levels at which $\tau$ chooses $P_0$),

($Z$ padded from the bottom),

(the levels that we have not reached).

Since the order of $\tau$ is $m-|Z| \ge t$, the set $V$ has exactly $m-t$ levels, namely, all the levels of $Q$ except the $t$ at which we saw a term other than $P_0$. We define $F_V$ by the same formula for every $V \subseteq Q$, with $F_V := Q$ if $V=Q$. It is the longest initial segment of $Q$ contained in $V$. ($Z$ stands for zero, since these are the levels where $\tau$ chooses $P_0$, and $F$ stands for free, since these are the levels where the box has a star.)

For $V \subseteq Q$ with $|V|=m-t$, let $\mathcal{B}_V$ be the set of the $9^t$ cubes $\pi$ with

Every cube in $\mathcal{B}_V$ is a box: it satisfies condition (i) of the definition because it has $|V|=m-t$ symbols $P_0$ or stars, and condition (ii) because its stars are below its symbols $P_0$, since $F_V$ is an initial segment of $Q$. When $V$ is the set defined above from a leaf $\tau$, the box of $\tau$ is the one in $\mathcal{B}_V$ with the terms of $\tau$ at the levels of $Q\setminus V$. Since a box of $\mathcal{B}_V$ has stars at the levels of $F_V$, it contains the leaves with any terms there, of which there are $10^{|F_V|}$. Notice that the set $\mathcal{B}_V$ depends on $w$ (since it depends on $Q$ and on the terms at the levels outside $Q$), but we omit this from the notation since $w$ will always be clear from context. We call the boxes in these sets, over all $V\subseteq Q$ with $|V|=m-t$, the boxes of $w$.

There is another way to describe the boxes of $w$, in terms of the leaves of order exactly $t$ contributing to $w$. For such a leaf, consider the lowest level of $Q$ at which it chooses a term other than $P_0$, and replace its $P_0$ by a star at every lower level of $Q$ (or at every level of $Q$, if $t = 0$). The result is a box of $w$, and each box of $w$ arises exactly once in this way. This means that $w$ has exactly $\alpha_t$ boxes.

Every leaf of order at least $t$ lies in exactly one box. A leaf contributing to $w$ that chooses $P_0$ at the set $Z$ of levels is a leaf of a box of $\mathcal{B}_V$ if and only if it chooses $P_0$ at every level of $V \setminus F_V$ and at no level of $Q \setminus V$, that is,

Furthermore, it is a leaf of exactly one of these boxes, namely, the one with its terms at the levels of $Q \setminus V$. We will see in the following lemma that, for $|Z| \leq m-t$, exactly one $V$ satisfies this condition, namely $Z$ padded from the bottom. Thus every leaf of order at least $t$ lands in exactly one box.

Lemma 27. Let $Q$ be a set of $m$ levels, and let $0 \leq t \leq m$. For every $Z \subseteq Q$ with $|Z| \leq m-t$, there is exactly one set $V \subseteq Q$ with $|V| = m-t$ and

namely $Z$ together with the $m-t-|Z|$ lowest levels of $Q \setminus Z$.

Proof. Let $V \subseteq Q$ with $|V| = m-t$ and $Z \subseteq V$. We show that $V \setminus F_V \subseteq Z$ holds if and only if $V \setminus Z$ consists of the $|V \setminus Z| = m-t-|Z|$ lowest levels of $Q \setminus Z$, which determines $V$. If $V = Q$, then both conditions hold: $F_V = Q$, so $V \setminus F_V = \varnothing\subseteq Z$, and $V \setminus Z$ is all of $Q \setminus Z$. Otherwise, let $\ell$ be the lowest level of $Q \setminus V$. Then $F_V$ consists of the levels of $Q$ below $\ell$, so $V \setminus F_V \subseteq Z$ says that every level of $V \setminus Z$ is below $\ell$. Since $Q \setminus Z$ is the disjoint union of $V \setminus Z$ and $Q \setminus V$, and every level of $Q \setminus V$ is at least $\ell$, this holds if and only if $V \setminus Z$ consists of the lowest levels of $Q \setminus Z$.

Lemma 28. Let $w$ be an output string, with inner set $Q$. The leaves of order at least $t$ contributing to $w$ are exactly the leaves of the boxes in $\bigcup_{V \subseteq Q,\ |V|=m-t} \mathcal{B}_V$, the union over all $V \subseteq Q$ with $|V|=m-t$, and each of them is a leaf of exactly one of these boxes. Hence

a sum of $\sum_{d=0}^{t-1}\alpha_d$ products and $\alpha_t$ values of boxes.

Proof. Every leaf of a box of $\mathcal{B}_V$ agrees with the private leaf of $w$ outside $Q$, so it contributes to $w$, and it chooses $P_0$ at most $|V|=m-t$ times, so its order is at least $t$. Conversely, a leaf $\tau$ contributing to $w$, of order at least $t$, chooses $P_0$ at a set $Z\subseteq Q$ of at most $m-t$ levels, and as we saw above, it is a leaf of a box of $\mathcal{B}_V$ if and only if $V\setminus F_V\subseteq Z\subseteq V$, and then of exactly one such box. By Lemma 27, exactly one $V$ satisfies this condition. The desired equation follows, since $(X_QY_Q)[w]$ is the sum of the products at the leaves contributing to $w$. Finally, $\alpha_d$ leaves of order $d$ contribute to $w$, and there are $\binom{m}{m-t}$ sets $V$ with $9^t$ boxes each, i.e., $\binom{m}{m-t}9^t=\binom{m}{t}9^t=\alpha_t$ boxes in all (as many as the leaves of order exactly $t$ contributing to $w$, although they cover its leaves of every order from $t$ to $m$).

The data structure, in terms of the parameters

In this subsection, we assemble the boxes into our data structure and bound its costs in terms of the parameters $m$, $L$, and $t$. We will then choose settings for the parameters in Section 4.4. We begin by recalling two notions from Section 2 that the data structure builds on.

The tiling of Section 2.3.4 cuts the product $XY$ into tiles. A row block is $N_0$ consecutive rows of $X$, a column block is $N_0$ consecutive columns of $Y$, and a tile is a band of $K_0 = \lfloor\sqrt{K}\rfloor$ consecutive row blocks together with a band of $K_0$ consecutive column blocks. One run of the recursion computes all $K_0^2$ block products of a tile at once, each as the product $X_QY_Q$ for its own subset $Q$. The encodings of Section 2.4.1 are the arrays of the $10^L$ numbers $\Phi_\tau(a)$, one for each leaf $\tau$, computed from the input array $a$ of a row band, and similarly the numbers $\Psi_\tau(b)$ computed from the input array $b$ of a column band. We compute an encoding once per band and share it among all the tiles that use that band, and the product at a leaf of a tile is the product of its two encoded numbers.

Our data structure consists of three components:

• the list of the $K_0^2$ subsets $Q$ assigned to the block products of a tile (recall that this list is the same for every tile),

• the encodings of all row bands and all column bands, and

• for every tile, the values of all its boxes (the boxes are the same strings in every tile, but their values depend on the input arrays of the tile).

To answer a query for an entry $(XY)[I,J]$, the data structure performs the following three steps:

1. Find the tile that contains $(I,J)$, the subset $Q$ of its block product, and the output string $w$ of $(I,J)$.

2. Add up the products at the leaves of order below $t$ contributing to $w$, by reading the two factors of each product from the encodings of the tile.

3. Add to this the values of the $\alpha_t$ boxes of $w$, which are stored for the tile.

By Lemma 28, the result is $(X_QY_Q)[w] = (XY)[I,J]$. We will give the details of the query, as well as the preprocessing, in the proof of Theorem 30 below.

The encodings are indexed by the leaves (Section 2.4.1), and for each tile we store its boxes, with their values, in a standard trie on their strings of $L$ symbols. Thus reading the product at a leaf, or looking up or inserting a box, takes $O(L)$ operations.

We will next show that there are at most $m + 1$ times as many boxes as there are leaves of order at least $t$, and that the dynamic program described in the technique overview computes their values (a box with $e$ stars being the sum of ten boxes with $e - 1$ stars).

Lemma 29. There are at most $(m + 1)\sum_{d=t}^{m}\beta_d$ boxes. Given the two encodings of a tile, we can compute the values of all these boxes, and store them in the trie for that tile, in $O(L)$ time and space per box.

Proof. The count. A box in which $f$ of the symbols are $P_0$ or stars is obtained from a leaf of order $m-f$ (namely the leaf we get by replacing its stars with $P_0$) by turning the $e$ lowest symbols $P_0$ of that leaf into stars, for some $e \leq f$. Thus, such a box can be specified by choosing the $f$ levels of these symbols, the number $e$, and a term other than $P_0$ at each of the other $L-f$ levels. There are hence $(f+1)\binom{L}{f}9^{L-f} = (f+1)\beta_{m-f}$ such boxes, and summing over $f \leq m-t$ gives the upper bound on the number of boxes.

The values. We compute the values of the boxes in increasing order of their number of stars. The boxes without stars are the leaves with at most $m-t$ symbols $P_0$, and the value of each of them is the product of its two numbers in the encodings.

For a box $\pi$ with $e \ge1$ stars, let $\ell$ be the highest level at which $\pi$ has a star. The leaves of $\pi$ are the leaves of the ten strings $\pi[\ell\leftarrow\lambda]$ obtained by replacing that star by a term $\lambda$, and each of these strings is again a box, with $e-1$ stars. Indeed, if $\lambda= P_0$, then the number of symbols that are $P_0$ or stars stays the same, and the new $P_0$ at level $\ell$ is above all the remaining stars. If $\lambda$ is any other term, then that number decreases by one, and every star is still below every $P_0$. Thus conditions (i) and (ii) in the definition of a box hold in both cases. (We expand the highest star because replacing a lower star by $P_0$ would leave a star above a $P_0$, which is not a box.) Hence we compute

with ten lookups in the trie, in $O(L)$ operations. Generating the boxes with $e$ stars and inserting them into the trie also takes $O(L)$ operations per box, and adds at most $L$ vertices per box.

See Figure 7 for an illustration of these counts (in that figure, the switching order $t$ is set to $m/9$). At every order $d$, the pruned recursion of Section 2.4 visits at most $\min\{|U|\alpha_d,\beta_d\}$ leaves of order $d$ in a tile whose wanted output strings form the set $U$ (see the proof of Lemma 11), the smaller of the two bounds at that order. The data structure cannot do this, since $U$ is not known in advance. Instead, a query computes $\alpha_d$ products at every order $d < t$, and looks up $\alpha_t$ boxes at the order $t$. The data structure also stores at most $(m+1)\sum_{d \ge t}\beta_d$ boxes per tile, even if they are not used by a query. The preprocessing time is therefore dominated by the sum $\sum_{d \ge t}\beta_d$, and hence by how fast the $\beta_d$ decay. Let

which is less than 1 since $L \ge10m$. It plays the role of $\frac{1}{2}$ in (5): for $1 \le d \le m$,

The closer the ratio $c=L/m$ is to 10, the closer $\rho$ is to 1, and the more boxes the preprocessing has to compute. This is one side of the trade-off mentioned at the beginning of this section between the time savings and how thin the product must be.

We can now bound the costs. The theorem below assumes that $N \ge\sqrt{K}N_0$. Since a tile spans $K_0N_0 \le\sqrt{K}N_0$ rows and as many columns of the $N \times N$ product $XY$, this ensures that a tile fits inside $XY$, so padding $N$ to a multiple of $K_0N_0$ at most doubles it. We checked the same requirement in Section 2.4.4, where it was the first of the two consequences of $N \ge D^{18}$. Here it will also follow, in Section 4.4, from the requirement that computing the encodings takes at most $N^2/D^\gamma$ time, since $10^L \ge M$.

Theorem 30. Let $m \ge1$, $D=4^m$, $L \ge10m$, $0 \le t \le m$, and $N \ge\sqrt{K}N_0$. Given as input matrices $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$, whose entries are integers of absolute value at most $N^{O(1)}$, we can preprocess them deterministically in time and space

After this, any single entry $(XY)[I,J]$ can be computed deterministically in $O\left(L\sum_{d=0}^{t}\alpha_d\right)$ time. In particular, for every set $W$ of positions of an $N \times N$ matrix, the entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in time

The constants hidden in the $O(\cdot)$ depend only on the exponent in $N^{O(1)}$.

Proof. Preprocessing. We tile the product as in Section 2.3.4, padding $N$ to a multiple of $K_0N_0$, which at most doubles it since $N\ge\sqrt{K}N_0\ge K_0N_0$, and we store the $K_0^2$ subsets, in $O(KL)\le O(10^L)$ operations. We then form and encode the input arrays of all row bands and column bands, at most $4N/(\sqrt{K}N_0)$ of them since $K_0\ge\sqrt{K}/2$. Each has at most $KN_0D\le7^L$ nonzero entries, so we form it in $O(L\cdot7^L)$ operations, and its encoding takes $O(10^L)$ operations (Section 2.4.1). This is the last term of (8). Finally, we compute the values of all the boxes of each of the at most $4N^2/M$ tiles, where $M=KN_0^2$ is the number of output entries of a tile. By Lemma 29, this takes $O(L(m+1)\sum_{d\ge t}\beta_d)$ time and space per tile, and $\sum_{d\ge t}\beta_d\le M\rho^t/(1-\rho)$ by (7), which gives the first term.

Query. Given $(I,J)$, we find its tile from the bands of row $I$ and column $J$, the subset $Q$ of its block product, and its output string $w$ (whose variables at the levels outside $Q$ we read off the row and column of $(I,J)$ within the block product), all in $O(L)$ operations (Section 2.4.4). We then compute $(X_QY_Q)[w]$ as the sum in Lemma 28, which has two parts.

The first part is the sum of the products at the leaves of order below $t$ contributing to $w$. Each of these leaves is obtained from the private leaf of $w$, which has $P_0$ at every level of $Q$, by picking fewer than $t$ of those levels and replacing $P_0$ with one of the nine other terms at each of them. We enumerate them, and for each such leaf $\tau$ we look up its two numbers $\Phi_\tau(a)$ and $\Psi_\tau(b)$ in the encodings of the tile and multiply them. The second part is the sum of the values of the boxes of $w$: for every $V\subseteq Q$ with $|V|=m-t$ and every box of $B_V$, we look up its value in the trie of the tile. In all, these are $\sum_{d=0}^{t-1}\alpha_d+\alpha_t$ numbers, and we find each of them in $O(L)$ operations, including forming its leaf or its box from the private leaf. Enumerating them also takes $O(L)$ operations per number. By Section 2.4.4, the sum is $(XY)[I,J]$. For a set $W$, we ask $|W|$ queries, which gives (9).

Word size. As in Section 2.4.4, every number in an encoding is a sum of entries of the input array with coefficients $0,\pm1$, and every value we compute is a sum of at most $10^m$ products of two such numbers. Since $N\ge N_0=3^{L-m}$ and $L\ge10m$, we have $10^L\le N^{5/2}$, so, as in Section 2.4.4, all the values are $N^{O(1)}$ in absolute value, and every integer we use has $O(\log N)$ bits.

Choosing the parameters

We now choose the parameters in Theorem 30 and prove the results of Section 4.1. In other words, we do the arithmetic and bookkeeping to derive the explicit exponents from Theorem 30. We first do the calculation in full for a single choice of the parameters, $L=21m$ and $t=\lceil m/9\rceil$, which proves Corollary 26. We then say what changes for other choices, which gives Theorems 24 and 25 and the entries of Table 2.

Proof of Corollary 26. Setting up. Let $m:=\lceil\log_4D\rceil$, and pad the inner dimension to $4^m<4D$ with zero columns of $X$ and zero rows of $Y$. This changes no entry of $XY$, and it changes $D$ by a factor less than 4, which only affects the constants. So from now on $D=4^m$, and since the original $D$ was larger than $4^{m-1}$, the assumption $N \ge D^{18}$ now reads $N \ge4^{18(m-1)}$. Let $L := 21m$ and $t := \lceil m/9\rceil$. We will apply Theorem 30 with these parameters, and we will verify its condition $N \ge\sqrt{K}N_0$ below, in the step on the encodings.

Boxes. We have $\rho= 9m/(L-m+1)=9m/(20m+1)<9/20$, so $1/(1-\rho)<2$, and since $t \ge m/9$, also $\rho^t \le(9/20)^{m/9}=D^{-\gamma}$, where $\gamma:= \ln(20/9)/(9\ln4)=0.0640\ldots$. Hence the first term of (8) is $O(m^2N^2/D^\gamma)$.

Queries. A query reads $\sum_{d\le t}\alpha_d$ numbers, and we bound this sum as in the proof of Lemma 11. Since $9^d=72^d8^{-d}\le72^t8^{-d}$ for $d\le t$, and $t<m/9+1$,

where $q:=\ln(72\cdot(9/8)^9)/(9\ln4)=0.4277\ldots$. Hence a query takes $O(LD^q)=O(mD^q)$ time.

Encodings. It remains to show that the last term of (8) is at most $N^2/D^\gamma$, and that $N\ge\sqrt{K}N_0$. Both follow from

the first by rearranging, and the second because $10^L\ge M=KN_0^2$ and $D^\gamma\ge1$. As in Section 2.4.4, we compare $m$-th powers by comparing their bases. We have $D^\gamma=(20/9)^{m/9}$, $10^L=10^{21m}$, and $N_0=3^{20m}$, and the standard bound $\binom{n}{k}\ge\frac{1}{n+1}\frac{n^n}{k^k(n-k)^{n-k}}$ gives $K=\binom{21m}{m}\ge\frac{1}{21m+1}(21^{21}/20^{20})^m$. Hence the left-hand side of (10) is at most $\sqrt{21m+1}\cdot\Lambda^m$, where

We now use the assumption $N\ge D^{18}$. Its base $4^{18}=6.871\ldots\cdot10^{10}$ is larger than $\Lambda$ by a factor of more than 1.63. Since $N\ge4^{18(m-1)}$, the inequality (10) therefore holds whenever $1.63^m\ge4^{18}\sqrt{21m+1}$, which is the case for all $m\ge60$. (For $m<60$, $D$ is bounded by a constant, and the corollary holds trivially.)

Conclusion. For $m\ge60$, Theorem 30 thus applies. The preprocessing takes $O(m^2N^2/D^\gamma)$ time and space, and a query takes $O(mD^q)$ time. Since $m=O(\log D)$ and the exponents $\gamma-0.063$ and $0.437-q$ are positive, we have $m^2=O(D^{\gamma-0.063})$ and $m=O(D^{0.437-q})$, which gives the bounds $O(N^2/D^{0.063})$ and $O(D^{0.437})$ of the corollary. The bound for a set $W$ follows by asking $|W|$ queries.

Other choices of the parameters. The numbers 21 and $\frac{1}{9}$ did not play a special role in this proof, and the same calculation works with $L=\lceil cm\rceil$ and $t=\lceil\theta m\rceil$ for any constants $c>10$ and $0<\theta<0.9$. This gives Corollary 31 below, which we use for the left half of Table 2, with the $c$ of the row and the $\theta$ that gives the $q$ of the column ($\theta=\frac{1}{9}$ for $q=0.43). To state it, let $H(x):=-x\ln x-(1-x)\ln(1-x)$ be the entropy function (with natural logarithms), let

be the general form of the base $\Lambda$ above, and let

Corollary 31. Let $c > 10$ and $0 < \theta< 0.9$, let $\rho_c := 9/(c - 1)$, and let

Given as input matrices $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$ whose entries are integers of absolute value at most $N^{O(1)}$, where $2 \le D \le N^\varepsilon$ and $\varepsilon< R_c(\gamma)$, with $R_c$ as in (11), we can preprocess them deterministically in $O(N^2 \log^2 D/D^\gamma)$ time and space, after which any single entry $(XY)[I,J]$ can be computed deterministically in $O(D^q \log D)$ time. The constants hidden in the $O(\cdot)$ depend on $c$, $\theta$, and $\varepsilon$.

Moreover, $R_c(0)$ increases to $\varepsilon^* = \ln4/(5 \ln10) = 0.1204\ldots$ as $c$ decreases to 10.

Proof. We repeat the proof of Corollary 26 with $L := \lceil cm\rceil$ and $t := \lceil\theta m\rceil$. After the padding, the assumption $D \le N^\varepsilon$ now reads $N \ge4^{(m-1)/\varepsilon}$. For the boxes, $\rho< \rho_c$ gives $\rho^t \le\rho_c^{\theta m} = D^{-\gamma}$. For the queries, the same calculation with $x := (1-\theta)/\theta$ in place of 8, where $9x > 1$ because $\theta< 0.9$, gives $\sum_{d \le t} \alpha_d \le(9x)^t(1+1/x)^m = O(D^q)$. For the encodings, the same bound on the binomial coefficient, $\binom{n}{k} \ge e^{nH(k/n)}/(n+1)$, shows that the left-hand side of (10) is $O(\sqrt{m}\cdot\Lambda^m)$. As $\varepsilon< R_c(\gamma)$ says that $\Lambda< 4^{1/\varepsilon}$, the inequality (10) thus holds once $m$ exceeds a constant depending on $c$, $\theta$, and $\varepsilon$, and for smaller $m$ the corollary again holds trivially.

Finally, a direct calculation shows that for $\gamma= 0$, $\ln\Lambda$ increases with $c$ and tends to $5\ln10$ as $c$ decreases to 10, by the identity $H(1/10) + (9/10)\ln9 = \ln10$. Hence $R_c(0)$ increases to $\ln4/(5\ln10) = \varepsilon^*$.

Where $\varepsilon^*$ comes from. The identity $H(1/10) + (9/10)\ln9 = \ln10$ says that $\binom{10m}{m}9^{9m}$, the largest term of the binomial expansion of $(1+9)^{10m}$, is $10^{10m}$ up to polynomial factors. So at $c=10$ a tile has $M \approx10^L$ output entries, as many as the recursion has leaves, and the encodings cost $N \cdot10^L/\sqrt{M} \approx N \cdot10^{L/2}$, which is below $N^2$ exactly when $N \ge10^{L/2} = D^{1/\varepsilon^*}$.

Proof of Theorem 24. Let $\varepsilon< \varepsilon^*$ and $q > 0$. Choose $c > 10$ with $R_c(0) > \varepsilon$ (Corollary 31), and then $\theta\in(0,0.9)$ small enough that $(H(\theta) + \theta\ln9)/\ln4 < q$ and $R_c(\theta\ln(1/\rho_c)/\ln4) > \varepsilon$. This is possible since, as $\theta\to0$, the first expression tends to 0 and the second to $R_c(0)$. Corollary 31 with these $c$ and $\theta$ gives the theorem.

Proof of Theorem 25. Apply Theorem 24 with $q := \kappa/2$, and ask one query for each position of $W$. As $|W| \le N^2/D^\kappa$, this takes $O(N^2 \log^2 D/D^\gamma+ |W|D^{\kappa/2}\log D) = O(N^2 \log^2 D/D^{\gamma'})$ time, where $\gamma' := \min\{\gamma,\kappa/2\}$.

With Corollary 31 in place of Theorem 24, the same argument gives explicit exponents, which we use for the right half of Table 2, with the $\theta$ at which $\gamma= \kappa- q$. The $\varepsilon$ of a row of the table is the smallest $R_c(\gamma)$ over its entries.

Corollary 32. Let $c$, $\theta$, $\gamma$, $q$, $X$, $Y$, $D$, and $\varepsilon$ be as in Corollary 31, and let $\kappa> 0$. For every set $W$ of at most $N^2/D^\kappa$ positions of an $N \times N$ matrix, the entries $(XY)[I,J]$, $(I,J) \in W$, can be computed deterministically in time

a saving of $D^{\min\{\gamma,\kappa-q\}}$, up to logarithmic factors, over the size $N^2$ of the product.

Further reductions

This section presents the remaining algorithmic results of the introduction: directed unweighted APSP (Section 5.1), the real-valued versions of 3SUM, APSP, and Exact Triangle (Section 5.2), the weighted $k$-Clique problems (Section 5.3), and the refutation of the rectangular hinted OMv conjectures of van den Brand, Nanongkai, and Saranurak [155] (Section 5.4). Each result but the last composes a known reduction with the algorithms of Sections 3 and 4; the last applies the data structure of Section 4 directly. We cite the reduction we use in each case and modify it only where our parameter regime requires it, saying what changes.

What we use. We recall the three results of Sections 3 and 4 that the reductions plug into.

• Exact Triangle (Section 3.2) asks whether a complete tripartite graph with parts $A,B,C$ of $n$ vertices and integer edge weights of absolute value at most $n^\nu$ has a triangle of weight zero. Theorem 19 solves it deterministically in $O(n^{3-\varepsilon_T}\log n) \le O(n^{3-\varepsilon_T})$ time, where $\varepsilon' = 0.00175$ and $\varepsilon_T = 0.0017$, and Theorem 22 derives from this a deterministic $\widetilde{O}(n^{3-\varepsilon'/3})$-time algorithm for the $(\min,+)$-product of two $n \times n$ matrices with integer entries of absolute value $n^{O(1)}$.

• Lop-AE-SparseTri$(n,D)$ and #Lop-AE-SparseTri$(n,D)$ (Definitions 13 and 14) are given a tripartite graph with parts $A$ and $B$ of $n$ vertices and a middle part of at most $D$ vertices, together with a set $W$ of query pairs $(a,b) \in A \times B$, and ask for every query pair in $W$ whether $a$ and $b$ have a common neighbor in the middle part, or for the counting version, how many they have. Corollary 16 solves both deterministically in $O(|W|D^{0.437} + n^2/D^{0.063})$ time when $n \ge D^{18}$. Theorem 17 reduces Exact Triangle to Lop-AE-SparseTri deterministically by hashing the weights modulo a prime.

• Corollary 26 is the data structure behind Corollary 16. Given $X \in\mathbb{Z}^{N \times D}$ and $Y \in\mathbb{Z}^{D \times N}$ with $N \ge D^{18}$ and entries of absolute value $N^{O(1)}$, it preprocesses them deterministically in $O(N^2/D^{0.063})$ time, after which it computes any entry of $XY$ in $O(D^{0.437})$ time. Section 5.2 applies it to matrices whose entries are not 0/1, and Section 5.4 uses its queries online.

The algorithms of Sections 5.1, 5.3, and 5.4 are deterministic. Those of Section 5.2 are randomized (Las Vegas), because the reductions of [?] that they use as black boxes are randomized.

Directed APSP with small integer weights

Let $\omega(1,r,1)$ denote the exponent of multiplying an $n \times n^r$ matrix by an $n^r \times n$ matrix, and let $\mu$ be the solution of $\omega(1,\mu,1)=1+2\mu$. Then $\mu\ge\frac{1}{2}$, since $\omega(1,\mu,1) \ge2$, and $\mu< 0.5275$ by the bounds of [10]. Zwick’s algorithm [172] solves APSP in directed unweighted graphs in $\widetilde{O}(n^{2+\mu})$ time, and the Directed Unweighted APSP hypothesis [66, ?, 82] asserts that no algorithm is polynomially faster. We refute this in Theorem 34.

We revisit Theorem 3.1 of Chan, Vassilevska Williams, and Xu [66], which is just a presentation of Zwick’s algorithm [172] as a reduction. Zwick’s algorithm is typically used in its randomized form because the proof is simple, so the reduction in [66] was stated in its randomized form. However, following an argument by Yuster and Zwick [171], Section 8, the reduction has a deterministic version as well:

Theorem 33 (Zwick’s algorithm as a deterministic reduction [172, 67]). If the $(\min,+)$-product of an $n \times n^\mu$ matrix by an $n^\mu\times n$ matrix, both with integer entries of absolute value at most $\widetilde{O}(n^{1-\mu})$, can be computed in $O(n^{2+\mu-\varepsilon_1})$ time for a constant $\varepsilon_1 > 0$, then APSP in directed $n$-vertex graphs with integer weights of absolute value polylog $n$ and no negative cycles can be solved in $O(n^{2+\mu-\varepsilon''})$ time for a constant $\varepsilon'' > 0$ that depends on $\varepsilon_1$ and $\mu$. If the algorithm for the product is deterministic, then so is the algorithm for APSP.

The randomized proof of [67] follows the randomized version of Zwick’s algorithm, which takes random samples of vertices to hit the long shortest paths, so the algorithm that it produces is randomized. Zwick [172], Sections 6 and 7 also gives a deterministic version of his algorithm, which constructs the hitting sets (his bridging sets) from witnesses of the products by a greedy hitting-set computation, within the same running time up to polylogarithmic factors, and Yuster and Zwick [170], Section 8 announce, without details, a faster deterministic construction of bridging sets for their distance oracle. Zwick’s construction lists a path with $s$ vertices for each of the $n^2$ pairs, which is more than the reduction can afford when $s$ is large. We therefore give a deterministic version of the reduction, which lists only walks with an endpoint in the previous hitting set, in the spirit of [170]. We believe that this is likely what Yuster and Zwick [170] had in mind.

Proof of Theorem 33. The deterministic algorithm follows the randomized reduction, with two changes. It constructs the hitting sets deterministically, by a hierarchy of bridging sets in the manner of Zwick [172], Section 6, and it computes the distances that this construction needs only to and from small sets of vertices, like the distance oracle of Yuster and Zwick [170]. The hypothesized algorithm for the $(\min,+)$-product is used only for the values of its products. (This is convenient but not essential: the witness-finding method of Seidel [145] and Alon, Galil, Margalit, and Naor [15, 91], derandomized by Alon and Naor [24] and applied to $(\min,+)$-products by Zwick [172], Lemma 3.3, uses the product algorithm only on submatrices, so any deterministic algorithm for the product can be made to return witnesses at a polylogarithmic cost.) The price is a factor $n^\gamma$ in the running time, for an arbitrarily small constant $\gamma> 0$, which the polynomial saving of the hypothesis absorbs. Throughout, $M$ is the largest absolute weight, which is polylogarithmic, and $\widetilde{O}$ hides polylogarithmic factors.

Bridges. For vertices $u$ and $v$, let $\delta(u,v)$ be the distance from $u$ to $v$, and let $\eta(u,v)$ be the smallest number of edges on a shortest path from $u$ to $v$. Every subpath of a shortest path with the fewest edges is again a shortest path with the fewest edges. For a threshold $t$ and a constant $C$, a set $B \subseteq V$ is a $(t,C)$-bridge if every pair $(u,v)$ with $\eta(u,v) \ge t$ has a shortest walk from $u$ to $v$ that passes through a vertex of $B$ and has at most $C\eta(u,v)$ edges. This is Zwick’s strong bridging set with a constant-factor slack in the number of edges, and with walks in place of paths because a shortest walk may run around zero-weight cycles. Let $R := \lceil n^\gamma\rceil$ (so $R$ is tiny!) and $t_i := \min\{R^i,n\}$ for $i = 0,1,\ldots,\lceil1/\gamma\rceil$. We shall construct $(t_i,C_i)$-bridges $B_i$ with $|B_i| = \widetilde{O}(n/t_i)$, starting with $B_0 = V$ and $C_0 = 1$, where $C_i$ is a constant that depends only on $i$. The slack $C_i$ grows by a constant factor from each level to the next. This is why the levels are a factor $n^\gamma$ apart, so that there are only $O(1/\gamma)$ of them, rather than a factor $3/2$ apart as in Zwick’s algorithm.

The tables $E_i(T,H)$. The distances to and from small sets of vertices are computed by the following routine. Given the bridges $B_0,\ldots,B_i$, a set $T$ of $\widetilde{O}(n/t_i)$ target vertices, and a horizon $H$ with $t_i \le H \le O(Rt_i)$, the routine returns two tables, $E_i(T,H)[V,T]$ and $E_i(T,H)[T,V]$, with three properties: every finite entry is the weight of a walk that the routine represents explicitly, as described next; every represented walk has $O(H)$ edges, with a constant that depends on $i$; and the entry for $(u,v)$ is $\delta(u,v)$ whenever $\eta(u,v) \le H$. The products in the routine are computed with the deterministic $(\min,+)$-product algorithm with witnesses of [172], which finds the witnesses deterministically by the method of [24]. A witness of an entry $(u,v)$ of a product is an index $k$ at which its minimum is attained, and we store it with the entry. The walk that the entry represents is then the concatenation of the walks represented by the entries $(u,k)$ and $(k,v)$ of the two factors, and it is listed, when needed, by following the stored witnesses recursively down to the entries of the weight matrix, which represent single edges; this takes $O(H)$ time, because the walk is a concatenation of $O(H)$ such entries. For a product in which one dimension is $n$ and the other two are at most $q$, with entries of absolute value at most $W$, that algorithm takes $\widetilde{O}(Wnq^{\omega-1})$ time, by cutting the product into $n/q$ square products.

The routine for level $0$ takes the weight matrix with zeroes on the diagonal and squares it $O(\log H)$ times, with the entries larger than $2HM$ replaced by $\infty$, and keeps the rows and columns of $T$. The routine for level $i \ge1$ first calls the routine for level $i-1$ with the target set $T \cup B_i$ and the horizon $L := 2C_i t_i$. This call is okay, because $|T \cup B_i| = \widetilde{O}(n/t_i) \le\widetilde{O}(n/t_{i-1})$ and $L = O(Rt_{i-1})$. Let $F$ be the pair of tables that it returns; they have exact distances for the pairs with $\eta\le L$, and they have the entries of $V \times(T \cup B_i)$ and $(T \cup B_i) \times V$. Let $Q$ be the smallest power of two above $H/t_i$, and set the diagonal of $F[B_i,B_i]$ to zero. The routine computes, by repeated squaring,

where $\star$ is the $(\min,+)$-product and the minimum is entrywise, and it computes $E_i(T,H)[T,V]$ symmetrically.

The entry for a pair $(u,v)$ with $v \in T$ and $\eta(u,v) \le H$ is exact. To see this, take a shortest path from $u$ to $v$ with the fewest edges, and cut it into blocks of $t_i$ edges followed by a remainder of fewer than $t_i$ edges. Each block is a shortest path with the fewest edges between its endpoints, so the bridge property of $B_i$ replaces it by a shortest walk between the same endpoints that has at most $C_i t_i$ edges and passes through a vertex of $B_i$; choose one such vertex in each block. The result is a shortest walk from $u$ to $v$ through at most $H/t_i$ chosen vertices of $B_i$, in which $u$, the chosen vertices, and $v$ are consecutively at most $L = 2C_i t_i$ edges apart. The table $F$ is exact on these consecutive pairs, so the product above finds $\delta(u,v)$. The represented walks have $O(H)$ edges, since those of $F$ have $O(L)$ edges and at most $Q+2 = O(H/t_i)$ of them are concatenated. Every product of the routine has one dimension $n$, the two others $\widetilde{O}(n/t_i)$, and entries $O(HM) = \widetilde{O}(Rt_i)$, so it takes $\widetilde{O}(Rt_i \cdot n \cdot(n/t_i)^{\omega-1}) \le\widetilde{O}(n^{\omega+\gamma})$ time, and the recursion has $O(1/\gamma)$ levels.

The stages. The algorithm maintains a matrix of distance estimates, initially the weight matrix with zeros on the diagonal, and runs one stage for each level $i = 0,1,\ldots$ until $t_{i+1} = n$. The stage for level $i$ computes the tables $E_i := E_i(B_i,C_i t_{i+1})$, then the $(\min,+)$-product

and takes its entrywise minimum with the matrix of estimates. Then, if $t_{i+1} < n$, it constructs $B_{i+1}$ from the walks represented in $E_i$, as described below. The product fixes every pair $(u,v)$ with $t_i \le\eta(u,v) \le t_{i+1}$: by the bridge property of $B_i$, gives $b \in B_i$ with $\delta(u,v) = \delta(u,b) + \delta(b,v)$ and $\eta(u,b),\eta(b,v) \le C_i\eta(u,v) \le C_i t_{i+1}$, so both $E_i[u,b]$ and $E_i[b,v]$ are exact. The stages thus find all distances, and every estimate is the weight of a walk, so none is too small.

The cost of the stages. Let $s := t_{i+1}$. The product of the stage for level $i$ has middle dimension $|B_i| = \widetilde{O}(n/t_i) = \widetilde{O}(Rn/s)$ and entries $\widetilde{O}(s)$. Cutting $B_i$ into $\widetilde{O}(R)$ groups of $\lceil n/s\rceil$ vertices turns it into $\widetilde{O}(n^\gamma)$ products of an $n \times\lceil n/s\rceil$ matrix by an $\lceil n/s\rceil\times n$ matrix with entries $\widetilde{O}(s)$. Up to this grouping, these are the products of the randomized reduction, and the analysis of [66] applies to them. Fix a small constant $a > 0$. When $s \le n^{1-\mu-a}$, the middle dimension $n^r := \lceil n/s\rceil$ has $r \ge\mu+ a$, and fast rectangular matrix multiplication computes the product in $\widetilde{O}(s n^{\omega(1,r,1)}) = O(n^{2+\mu-ca})$ time [172], for a constant $c > 0$, by the convexity of $\omega(1,\cdot,1)$, since $\omega(1,\mu,1) = 1 + 2\mu$ and $\omega< 2 + \mu$. When $s \ge n^{1-\mu+a}$, the naive algorithm computes it in $O(n^2 \cdot n/s) = O(n^{2+\mu-a})$ time. In between, the middle dimension is at most $n^{\mu+a}$ and the entries are $\widetilde{O}(n^{1-\mu+a})$; we cut the middle dimension into at most $n^a$ pieces of $n^\mu$ columns and pad each piece to an instance with $n^{1+O(a)}$ rows, so that its entries fit the range of the hypothesis, and the hypothesized algorithm computes it in $O(n^{2+\mu-\varepsilon_1+O(a)})$ time. With $a$ small enough, the product of every stage takes $O(n^{2+\mu+\gamma-\varepsilon_2})$ time, for a constant $\varepsilon_2 > 0$ that depends on $\varepsilon_1$ and $\mu$.

The next bridge. Let $s := t_{i+1} < n$. The bridge $B_{i+1}$ is read off the walks represented in $E_i = E_i(B_i, C_i s)$, and this is where the small target set pays off. These are $2n|B_i|$ walks with $O(s)$ edges each, so listing them all takes $\widetilde{O}(n|B_i|s) = \widetilde{O}(n^2R) \le\widetilde{O}(n^{2+\gamma})$ time, whereas listing a path for every pair of vertices, as Zwick’s construction does, would take $O(n^2s)$ time, which exceeds our budget at the large scales. From each listed walk, erase the cycles, and keep the resulting path if it has at least $s/2$ edges. Let $B_{i+1}$ be a greedy hitting set of the vertex sets of the kept paths, each of which has more than $s/2$ vertices. It has $O((n/s)\log n)$ vertices by the analysis of the greedy heuristic [122, 56], and it is computed in time linear in the total size of the sets, up to a logarithmic factor.

This $B_{i+1}$ is an $(s,C_{i+1})$-bridge for a constant $C_{i+1}$ that depends only on $i$. Consider first a pair $(x,y)$ with $\eta(x,y) = s$. The bridge $B_i$ gives $b$ with $\delta(x,y) = \delta(x,b) + \delta(b,y)$ and $\eta(x,b) + \eta(b,y) \le C_i s$, so the entries of $E_i$ for $(x,b)$ and $(b,y)$ are exact, and the two walks that represent them are shortest walks. A cycle on a shortest walk has weight zero, so the two paths obtained by erasing the cycles are still shortest, and their concatenation is a shortest walk from $x$ to $y$ with $O(s)$ edges. This walk has at least $s$ edges, since otherwise erasing its cycles would give a shortest path from $x$ to $y$ with fewer than $\eta(x,y)$ edges. So one of the two paths has at least $s/2$ edges and was kept, and $B_{i+1}$ contains one of its vertices, which lies on the concatenated walk. For a pair with $\eta(x,y) > s$, apply this to the first $s$ edges of a shortest path from $x$ to $y$ with the fewest edges, and keep the rest of that path; the walk obtained has at most $C_{i+1}s + \eta(x,y) - s \le C_{i+1}\eta(x,y)$ edges.

Running time. There are $O(1/\gamma)$ stages. Each computes its tables in $\widetilde{O}(n^{\omega+\gamma})$ time, its product in $O(n^{2+\mu+\gamma-\varepsilon_2})$ time, and its bridge in $\widetilde{O}(n^{2+\gamma})$ time. Choosing $\gamma$ smaller than both $\varepsilon_2/2$ and $(2+\mu-\omega)/2$ gives the running time $O(n^{2+\mu-\varepsilon''})$ for a constant $\varepsilon'' > 0$, and every step is deterministic.

We compose the known reduction with the deterministic $(\min,+)$-product algorithm of Theorem 22.

Theorem 34 (Directed unweighted APSP). There is a constant $\varepsilon'' > 0$ such that APSP in directed unweighted $n$-vertex graphs, and more generally in directed $n$-vertex graphs with integer weights of absolute value at most $(\log n)^c$ and no negative cycles, for any constant $c$, can be solved in $O(n^{2+\mu-\varepsilon''})$ time by a deterministic algorithm.

Proof. Cut the $n \times n^\mu$ matrix into $n^{1-\mu}$ blocks of $n^\mu$ consecutive rows, and the $n^\mu\times n$ matrix into $n^{1-\mu}$ blocks of $n^\mu$ consecutive columns. The $(\min,+)$-product then consists of $n^{2-2\mu}$ products of two $n^\mu\times n^\mu$ matrices, one for each pair of blocks. Their entries have absolute value $\widetilde{O}(n^{1-\mu}) \le\widetilde{O}(n^\mu)$, because $\mu\ge\frac{1}{2}$, so Theorem 22 computes each of them in $\widetilde{O}((n^\mu)^{3-\varepsilon'/3})$ time, and all of them in $\widetilde{O}(n^{2+\mu-\mu\varepsilon'/3})$ time. This is $O(n^{2+\mu-\varepsilon_1})$ for any constant $\varepsilon_1 < \mu\varepsilon'/3$, so [8] applies, in the deterministic form of Theorem 33, since the algorithm of Theorem 22 is deterministic.

3SUM, APSP, and Exact Triangle with real inputs

Theorem 35 (Real inputs). In the real RAM in which the only operations allowed on real numbers are comparisons, additions, and subtractions, Las Vegas algorithms solve the following problems on real inputs: 3SUM on $n$ numbers in $O(n^{1.998})$ expected time, and the $(\min,+)$-product of two $n \times n$ matrices, APSP with no negative cycles, and Exact Triangle in $O(n^{2.998})$ expected time. The bounds also hold with probability $1-n^{-c}$ for any constant $c$.

The reductions are by Chan, Vassilevska Williams, and Xu [?], who reduce the real-valued problems to the counting version of All-Edges Sparse Triangle on graphs with $m$ edges (and real APSP also to the Boolean version). Our version of All-Edges Sparse Triangle in Corollary 16 needs lopsided sparse graphs. To deal with this slight change, we follow the reductions of [?] down to the point where the All-Edges Sparse Triangle oracle is called, and change the call. At that point, the reductions only need the oracle to compute comparison counts, which we define next. These come from Fredman’s trick [87]: $a' + b' < a + b$ if and only if $a' - a < b - b'$. This turns a comparison of two sums into a comparison of a number that depends only on the row against a number that depends only on the column.

We first define comparison counts, then restate the reductions of [?] as reductions to comparison counts (Lemma 36), and then solve the comparison counts problem: following [?], Lemma 37 reduces it to the wanted entries of one thin matrix product using Matoušek’s technique for the dominance product [129], and Corollary 38 applies Corollary 26 to that product. The proof of Theorem 35 then adds up the running times. The only randomized steps are in the reductions of [?], and the algorithms never output a wrong answer, so they are Las Vegas.

Comparison counts. We are given $n$ row lists: for every $r \in[n]$, a list $\mathcal{R}_r$ of at most $d$ real numbers $x$, where each $x \in\mathcal{R}_r$ has a color $\chi(x) \in[d]$ associated with it.

Similarly, we are given $n$ column lists: for every $c \in[n]$, a list $\mathcal{C}_c$ of at most $d$ real numbers $y$, each $y \in\mathcal{C}_c$ again with a color $\chi(y) \in[d]$ associated with it.

Given a set $P \subseteq[n]^2$ of (row, column) pairs, we want, for every $(r,c) \in P$, the number

The reductions. We follow the reductions of [?] down to their two counting problems [?], Problems 4.1 and 4.5, which we define in the proof of Lemma 36 below, and turn every call of a counting problem into comparison counts, in place of the triangle-counting instances of [?], Lemmas 4.3 and 4.7. Throughout, $d \in[n]$ is a parameter, which the proof of Theorem 35 chooses.

Lemma 36 (after Chan, Vassilevska Williams, and Xu [?]). Let $1 \le d \le n$.

• (a) The $(\min,+)$-product of two $n \times n$ real matrices, APSP with real edge weights and no negative cycles on $n$ vertices, and Exact Triangle with real weights on $n$ vertices reduce, with Las Vegas randomization, to $\widetilde{O}(n/d)$ counting calls and $\widetilde{O}(n^3/d)$ further time. A counting call consists of $d$ comparison counts, each with $n$ row lists and $n$ column lists of at most $d$ numbers, with at most one number of each color in a list, and with pair sets $P_1,\ldots,P_d$ satisfying $\sum_k |P_k| = n^2$; forming its lists costs $O(nd^2)$ subtractions.

• (b) $3$SUM on $n$ real numbers reduces, with Las Vegas randomization, to polylogarithmically many comparison counts, each with $n$ row lists and $n$ column lists of at most $d$ numbers, all of the same color, and a pair set $P$ of size $O(n^2/d)$, and $\widetilde{O}(n^2/d)$ further time; forming the lists of a comparison count costs $O(nd)$ subtractions.

A comparison count may have $\le$ in place of $<$. The only operations on real numbers performed by the reductions are comparisons, additions, and subtractions.

Proof. (a) $(\min,+)$-product, APSP, and Exact Triangle. Let $A,B,C \in\mathbb{R}^{n \times n}$. As in [67], Section 3.2, cut $A$ into $n/d$ blocks $A' \in\mathbb{R}^{n \times d}$ of $d$ consecutive columns, and $B$ into the corresponding blocks $B' \in\mathbb{R}^{d \times n}$ of $d$ consecutive rows. For each pair $(A',B')$ we solve [67], Problem 3.2: find, for every $(i,j)$, the predecessor and the successor of $C[i,j]$ among the $d$ sums $A'[i,k] + B'[k,j]$, where a sum equal to $C[i,j]$ counts as its predecessor. For Exact Triangle, $A[i,k] := w(i,k)$, $B[k,j] := w(k,j)$, and $C[i,j] := -w(i,j)$, and there is a triangle of weight zero if and only if, in some block, some $C[i,j]$ is its own predecessor. For the $(\min,+)$-product of $A$ and $B$, we take every $C[i,j]$ smaller than all the sums, for instance $C[i,j] := \min_k A[i,k] + \min_k B[k,j] - 1$. Then the successor of $C[i,j]$ in a block is the smallest sum of the block, and the product is the entrywise minimum over the $n/d$ blocks, which costs $O(n^3/d)$ comparisons. APSP with real weights and no negative cycles can be computed by successive squaring of the weight matrix, performing $\lceil\log_2 n \rceil$ such products.

Fix a pair $(A',B')$. We solve the following counting problem [67], Problem 4.1: we are given a pivot $k_{ij} \in[d]$ for every $(i,j) \in[n]^2$ and a subset $S \subseteq[d]$, and we must count, for every $(i,j)$, the indices $k' \in S$ with $A'[i,k'] + B'[k',j] < A'[i,k_{ij}] + B'[k_{ij},j]$; a call of the problem may also have $\le$ in place of $<$.

By Lemmas 3.4 and 4.2 of [67], the problem for one pair $(A',B')$ reduces, with Las Vegas randomization, to polylogarithmically many calls of the counting problem, and $\widetilde{O}(n^2)$ further time: Lemma 3.4 finds the predecessors and successors by a randomized search that repeatedly asks, for every $(i,j)$, for a sum strictly between two given sums, and Lemma 4.2 finds these sums by counting, for each of the two given sums, the sums below it, and by sampling at random. Lemma 4.3 of [67] turns a call into $O(d)$ triangle-counting instances; we turn it into $d$ comparison counts instead.

Following [67], we group the pairs by their pivot, $P_k := \{(i,j): k_{ij} = k\}$. For $k \in[d]$, let the row list $\mathcal{R}_i$ consist of the numbers $A'[i,k'] - A'[i,k]$, $k' \in S$, each with color $k'$, and let the column list $\mathcal{C}_j$ consist of the numbers $B'[k,j] - B'[k',j]$, $k' \in S$, again with color $k'$.

By Fredman’s trick, $A'[i,k'] + B'[k',j] < A'[i,k] + B'[k,j]$ if and only if $A'[i,k'] - A'[i,k] < B'[k,j] - B'[k',j]$, where both numbers in the latter comparison have color $k'$.

Thus the count of a pair $(i,j) \in P_k$ is the comparison count $\gamma(i,j)$ of these lists with $P := P_k$. Forming the lists costs us $O(nd)$ subtractions for each $k$. Over the $n/d$ pairs $(A',B')$, and the $\lceil\log_2 n \rceil$ products in the case of APSP, this gives $\widetilde{O}(n/d)$ counting calls and $\widetilde{O}(n^3/d)$ additional time, which includes the $O(n^3/d)$ comparisons of the entrywise minimum.

(b) 3SUM. Let $A,B,C$ be sets of $n$ reals, the problem asks whether $a+b=c$ for some $a\in A$, $b\in B$, and $c\in C$. The 3SUM version with one set and $a+b+c=0$ reduces to this version in the standard way.

Now, sort $A$ and $B$, and cut them into $n/d$ consecutive blocks $A_1,\ldots,A_{n/d}$ and $B_1,\ldots,B_{n/d}$ of $d$ numbers each, writing $A_i[k]$ for the $k$-th number of $A_i$. Here the counting problem [67] is the following: we are given a set $Q$ of $O(n^2/d)$ quadruples $(i,j,k,\ell)$, and we must count, for every $(i,j,k,\ell)\in Q$, the pairs $(k',\ell')\in[d]^2$ with $A_i[k']+B_j[\ell']<A_i[k]+B_j[\ell]$; again, a call may have $\leq$ in place of $<$, and it may restrict $k'$ and $\ell'$ to subsets of $[d]$. By Section 3.3 and Lemmas 3.9 and 4.6 of [67], 3SUM reduces, with Las Vegas randomization, to polylogarithmically many calls of this counting problem, and $\widetilde{O}(n^2/d)$ additional time. By Fredman’s trick, the condition is $A_i[k']-A_i[k]<B_j[\ell]-B_j[\ell']$. So a call is one comparison count: the rows are the $n$ pairs $(i,k)$ and the columns are the $n$ pairs $(j,\ell)$, the row list $R_{(i,k)}$ consists of the $d$ numbers $A_i[k']-A_i[k]$, $k'\in[d]$, the column list $C_{(j,\ell)}$ of the $d$ numbers $B_j[\ell]-B_j[\ell']$, $\ell'\in[d]$, all numbers have the same color, and $P:=Q$. A restriction of $k'$ or of $\ell'$ deletes numbers from the lists. Forming the lists can be done using $O(nd)$ subtractions.

Solving comparison counts. In the comparison counts of Lemma 36(a), every list has at most one number of each color, so that $\gamma(r,c)$ is the number of colors $k$ with $A[r,k]<B[k,c]$, where $A[r,k]$ is the number of color $k$ in $R_r$ and $B[k,c]$ the one in $C_c$ (a color missing from either list is skipped); this is the dominance product of the $n\times d$ matrix $A$ and the $d\times n$ matrix $B$ at the positions in $P$. In those of Lemma 36(b), all numbers have the same color, and $\gamma(r,c)$ is the number of pairs $x<y$ between the two lists. Merging the two sorted lists of every pair in $P$ takes $O(|P|d)$ time. To improve on this, we utilize Matoušek’s technique for the dominance product [129], which counts the pairs that lie in different blocks of the sorted order by one thin matrix product, and the pairs in the same block directly.

Lemma 37 (after Matoušek [129]). Let $s\geq1$ be an integer and $D'' := \lceil2nd/s\rceil+d$. Given row lists $R_r$ and column lists $C_c$ as above and a set $P\subseteq[n]^2$, we can build matrices $X\in\{0,\ldots,d\}^{n\times D''}$ and $Y\in\{0,\ldots,d\}^{D''\times n}$ and compute integers $\gamma_2(r,c)$, $(r,c)\in P$, with

deterministically in $O(nD''+(|P|+nds+ds^2)\log n)$ time. The only operations on real numbers are the comparisons made when sorting the numbers of each color. The same holds with $\leq$ in place of $<$ in the definition of $\gamma$.

Proof. For each color $k$, sort the numbers of color $k$ of all the lists, and break ties so that, among equal numbers, those from column lists precede those from row lists. Let $\operatorname{pos}(x)$ be the position of $x$ in this order. For $x$ from a row list and $y$ from a column list of the same color, $x<y$ holds if and only if $\operatorname{pos}(x)<\operatorname{pos}(y)$. (For $\leq$, we reverse the tie-breaking rule.) Cut the order into blocks of $s$ consecutive positions, and let $\operatorname{blk}(x):=\lfloor\operatorname{pos}(x)/s\rfloor$. Then $x<y$ holds if and only if either $\operatorname{blk}(x)<\operatorname{blk}(y)$, or $\operatorname{blk}(x)=\operatorname{blk}(y)$ and $\operatorname{pos}(x)<\operatorname{pos}(y)$. As in Matoušek’s algorithm, we count the pairs of the first kind, which lie in different blocks, by a matrix product, and the pairs of the second kind, which lie in the same block, by brute force.

Different blocks. The lists contain at most $2nd$ numbers in all, so there are at most $2nd/s+d\leq D''$ pairs (color, block), and we index the columns of $X$ and the rows of $Y$ by these pairs, padding with zeros if needed. Let $X[r,(k,\beta)]$ be the number of numbers of color $k$ in $R_r$ that lie in block $\beta$, and let $Y[(k,\beta),c]$ be the number of numbers of color $k$ in $C_c$ that lie in a block after $\beta$. Then $(XY)[r,c] = \sum_{(k,\beta)} X[r,(k,\beta)]Y[(k,\beta),c]$ is the number of pairs $(x,y)$ of the same color with $\operatorname{blk}(x) < \operatorname{blk}(y)$, as required, and the entries of $X$ and $Y$ are between 0 and $d$ because a list has at most $d$ numbers. We use $O(nd\log n)$ comparisons to sort and then fill in $X$ and $Y$ in $O(nD^{\prime\prime})$ time.

Same block. For every color $k$ and block $\beta$, let $I$ be the numbers $x$ of row lists in the block, each with its row, and let $J$ be the numbers $y$ of column lists in the block, each with its column, so that $|I|+|J|\le s$. For every $(x,y)\in I\times J$ with $\operatorname{pos}(x)<\operatorname{pos}(y)$ whose row and column form a pair $(r,c)\in P$ (this can be tested via a binary search in $P$ after it is sorted), we add 1 to $\gamma_2(r,c)$. Then $\gamma_2$ counts the pairs of the second kind. Since $|I||J|\le s^2/4$ and there are at most $2nd/s+d$ blocks, we enumerate at most $(2nd/s+d)s^2/4\le nds+ds^2$ pairs $(x,y)$. Together with sorting $P$ and outputting the $|P|$ integers $\gamma_2(r,c)$, this takes $O((|P|+nds+ds^2)\log n)$ time.

Corollary 38. Let $d\le n^{1/40}$, and let $n$ be larger than a suitable constant. Given row lists $\mathcal{R}_r$ and column lists $\mathcal{C}_c$ as above and a set $P\subseteq[n]^2$, the counts $\gamma(r,c)$, $(r,c)\in P$, with $<$ or with $\le$, can be computed deterministically in

time, by an algorithm whose only operations on real numbers are comparisons.

Proof. Apply Lemma 37 with $s:=\lceil n^{1-1/440}/d\rceil$, so that $D^{\prime\prime}\le2d^2n^{1/440}+d+1\le3d^2n^{1/440}$, and pad the product to the inner dimension $D^*:=\lceil3d^2n^{1/440}\rceil$ with zero columns of $X$ and zero rows of $Y$. Since $d\le n^{1/40}$, we have $D^*\le4n^{23/440}$, so $(D^*)^{18}\le4^{18}n^{0.941}\le n$, and Corollary 26 with $N:=n$ and $D:=D^*$ computes $(XY)[r,c]$ for all $(r,c)\in P$ deterministically in $O(|P|(D^*)^{0.437}+n^2/(D^*)^{0.063})$ time. Here $(D^*)^{0.437}\le4^{0.437}n^{0.437\cdot23/440}=O(n^{0.0229})$ and $(D^*)^{0.063}\ge(3d^2)^{0.063}\ge d^{0.126}$.

Since $nD^*=O(n^{1.06})$, $nds\le2n^{2-1/440}$ and $ds^2\le4n^{2-1/440}$, the time of Lemma 37 itself is

Proof of Theorem 35. Let $d:=\lfloor n^{1/40}\rfloor$, apply Lemma 36 with this $d$, and answer every comparison count by Corollary 38.

$(\min,+)$-product, APSP, and Exact Triangle. By Corollary 38, the $d$ comparison counts of a counting call cost

since $\sum_k|P_k|=n^2$, $d\cdot n^2/d^{0.126}\le n^{2+0.874/40}\le n^{2.0219}$, and $d\cdot n^{2-1/440}\le n^{2+1/40-1/440}\le n^{2.0228}$; forming the lists costs $O(nd^2)=O(n^{1.05})$. So the $\widetilde{O}(n/d)$ counting calls of Lemma 36(a) cost

expected time, and its additional time is $\widetilde{O}(n^3/d)=\widetilde{O}(n^{2.975})$.

3SUM. By Corollary 38, a comparison count of Lemma 36(b) costs

so 3SUM can also be solved in $\widetilde{O}(n^{1.9979}) \le O(n^{1.998})$ expected time.

Las Vegas, and high probability. The reductions of [67] never return a wrong answer and Corollary 38 is deterministic, so the output is always correct, and only the running time is random. The bounds also hold with probability $1-n^{-c}$ for any constant $c$: we stop a run that exceeds twice the bound on its expected time and start again, and $c\log^2 n$ runs all fail with probability at most $n^{-c}$.

Zero-Weight, Min-Weight, and Max-Weight $k$-Clique

The classical reduction of Nešetřil and Poljak [132] from $k$-Clique to triangle detection turns the weighted $k$-Clique problems into Exact Triangle instances with $n^{\lfloor k/3\rfloor}$ vertices per part (after fixing $k-3\lfloor k/3\rfloor$ vertices when $3\nmid k$), so Theorem 19 gives:

Corollary 39 (Weighted $k$-Clique). Let $k \ge3$ and $\nu\ge1$ be constants. Given a complete $k$-partite graph with parts of $n$ vertices and integer edge weights of absolute value at most $n^\nu$, deterministic algorithms decide in $O(n^{k-\varepsilon_T\lfloor k/3\rfloor})$ time whether some $k$-clique, with one vertex in each part, has total edge weight zero, and find a $k$-clique of minimum, or of maximum, total edge weight.

Proof. We apply the classical reduction from $k$-Clique to triangle detection of Nešetřil and Poljak [132]. When $k$ is divisible by 3, the reduction is simple. Split the $k$ parts into three groups of $k/3$ parts each, and build a tripartite graph $H$ whose vertices in the $i$th part (for $i=1,2,3$) are the $k/3$-cliques of the original graph with one vertex in each part of the $i$th group, with an edge between two such cliques when together they form a $2k/3$-clique. Since the original graph is complete $k$-partite, every choice of one vertex from each part of a group is a $k/3$-clique, and any two such cliques from different groups form a $2k/3$-clique, so $H$ is simply the complete tripartite graph with $n^{k/3}$ vertices per part. The triangles of $H$ correspond bijectively to the $k$-cliques of the original graph with one vertex in each part.

When $k$ is not divisible by 3, let $t := k-3\lfloor k/3\rfloor\in\{1,2\}$. We enumerate the $n^t$ ways of fixing one vertex in each of the first $t$ parts, and for each of them we build $H$ as above from the remaining $3\lfloor k/3\rfloor$ parts; now the triangles of $H$ correspond to the $k$-cliques through the $t$ fixed vertices. In both cases $H$ has $N := n^{\lfloor k/3\rfloor}$ vertices per part, and we index its parts by $i \in\mathbb{Z}_3$.

We give $H$ edge weights so that a triangle and its $k$-clique have the same weight. Let $U$ be a vertex of $H$ in part $i$ and $V$ a vertex in part $i+1$. The weight of the edge $UV$ of $H$ is the total weight of the edges of the original graph between a vertex of $U$ and a vertex of $V$, of the edges inside $U$, of the edges between $U$ and the fixed vertices, and, if $i=1$ and $t=2$, of the edge between the two fixed vertices. The edges inside $V$ and the edges between $V$ and the fixed vertices are not counted here; they are counted in the edge from $V$ to part $i+2$. In this way every edge of a $k$-clique is counted exactly once over the three edges of its triangle, so there is no double-counting and the weights agree. The weight of an edge of $H$ has absolute value at most $\binom{k}{2}n^\nu\le N^{2\nu}$ for $n \ge k^2$.

Theorem 19 decides whether $H$ has a triangle of weight zero in $O(N^{3-\varepsilon'\log N})$ time, with the constants $\varepsilon' > \varepsilon_T$ of that theorem. By [158], Theorem 3.3, finding a maximum weight triangle costs $O(\log^2 n)$ times as much, and a minimum weight triangle is a maximum weight triangle for the negated weights. Since $\varepsilon' > \varepsilon_T$, the logarithmic factors are absorbed, so each graph $H$ costs $O(N^{3-\varepsilon_T})$ time, and the $n^t$ graphs together cost $O(n^tN^{3-\varepsilon_T}) = O(n^{k-\varepsilon_T\lfloor k/3\rfloor})$.

Consequently, the Zero-Weight, Min-Weight, and Max-Weight $k$-Clique hypotheses [21], [125, 36] with polynomially bounded weights are refuted, also on general graphs, by the standard reduction to the complete $k$-partite form.

Three conjectures of van den Brand, Nanongkai, and Saranurak

Van den Brand, Nanongkai, and Saranurak [155] propose three “hinted” variants of the Online Matrix–Vector conjecture [100], which they use to prove matching lower bounds for the running times of their dynamic matrix inverse algorithms. All three are over the Boolean semiring, allow polynomial preprocessing time in the first phase, and have parameters $t=n^{\tau}$ and $t_i=n^{\tau_i}$ for constants $\tau,\tau_i \in(0,1)$. We write $\omega(a,b,c)$ for the exponent of multiplying an $n^a \times n^b$ matrix by an $n^b \times n^c$ matrix. Each conjecture asserts that no algorithm simultaneously beats all of its listed bounds, for any $\varepsilon> 0$.

• $v$-hinted Mv (Definition 5.1 and Conjecture 5.2 of [155]). Phase 1: an $n \times t$ matrix $M$. Phase 2: a $t \times n$ matrix $V$. Phase 3: an index $i$, after which the algorithm outputs $MV_{[n],i}$, the product of $M$ and column $i$ of $V$. Bounds conjectured to be impossible to achieve simultaneously: $n^{\omega(1,1,\tau)-\varepsilon}$ for Phase 2 and $n^{1+\tau-\varepsilon}$ for Phase 3.

• $Mv$-hinted Mv (Definition 5.6 and Conjecture 5.7 of [155]). Phase 1: $N \in\{0,1\}^{n \times n}$ and $V \in\{0,1\}^{t \times n}$. Phase 2: a vector $I \in[n]^t$ of column indices. Phase 3: an index $j$, after which the algorithm outputs $N_{[n],I}V_{[t],j}$, where $N_{[n],I}$ is the $n \times t$ matrix whose $k$-th column is column $I_k$ of $N$. Bounds conjectured to be impossible to achieve simultaneously: $n^{\omega(1,\tau,1)-\varepsilon}$ for Phase 2 and $n^{1+\tau-\varepsilon}$ for Phase 3.

• $uMv$-hinted $uMv$ (Definition 5.11 and Conjecture 5.12 of [155]). Phase 1: $U \in\{0,1\}^{n \times t_1}$, $N \in\{0,1\}^{n \times n}$, and $V \in\{0,1\}^{t_2 \times n}$. Phase 2: $I \in[n]^{t_1}$. Phase 3: $J \in[n]^{t_2}$. Phase 4: indices $i$ and $j$, after which the algorithm outputs $(UN_{I,J}V)_{i,j}$, where $N_{I,J}$ is the $t_1 \times t_2$ submatrix of $N$ with the rows $I$ and the columns $J$. Bounds conjectured to be impossible to achieve simultaneously: $n^{\omega(1,\tau_1,1)-\varepsilon}$ for Phase 2, $n^{\omega(\tau_2,\tau_1,1)-\varepsilon}$ for Phase 3, and $n^{\tau_1+\tau_2-\varepsilon}$ for Phase 4.

The two obvious algorithms are to multiply everything out with fast matrix multiplication as soon as the hint arrives, or to wait and compute one matrix–vector product at the end, and the conjectures say that no algorithm is polynomially faster than both of these. The main lower bounds of [155], $n^{1.528}$ for dynamic matrix inverse with column updates and row queries and $n^{1.407}$ with element updates and element queries, use the conjectures at $\tau\approx0.53$ and at $\tau_1+\tau_2 \approx1.41$, where they match the algorithms of [155]. For thin hints, the conjectures are false, because the data structure of Corollary 26 is faster than both obvious algorithms.

Corollary 40 (Hinted OMv). Conjectures 5.2 and 5.7 of [155] fail for every $0 < \tau< 1/18$: Phase 2 takes $O(n^{2-0.063\tau})$ time and Phase 3 takes $O(n^{1+0.437\tau})$ time. Conjecture 5.12 fails for every $0 < \tau_1 < \tau_2/18$: Phase 3 takes $O(n^{1+\tau_2-0.063\tau_1})$ time and Phase 4 takes $O(n^{\tau_2+0.437\tau_1})$ *time, with polynomial Phase 1 and a Phase 2 that only stores $I$.

More generally, Conjectures 5.2 and 5.7 fail for every $0 < \tau< \varepsilon^* = 0.1204\ldots$, and Conjecture 5.12 fails for every $0 < \tau_1 < \varepsilon^*\tau_2$. In each case the phases beat the conjectured bounds whatever the value of $\omega$, since $\omega(1,1,\tau)=\omega(1,\tau,1)\ge2$ and $\omega(\tau_2,\tau_1,1)\ge1+\tau_2$ because of the input and output sizes.

Proof. We compute over $\mathbb{Z}$ with $0/1$ matrices; a Boolean entry is 1 exactly when the corresponding integer entry is positive. Let $D := t = n^{\tau}$ be the inner dimension. For $\tau< 1/18$ we have $D^{18} = n^{18\tau} \le n$, so Corollary 26 applies with $N = n$.

Conjecture 5.2 of [155]. In Phase 2, we preprocess $X := M$ and $Y := V$ by Corollary 26, in $O(n^2/D^{0.063}) = O(n^{2-0.063\tau})$ time. In Phase 3, the $n$ entries of column $i$ of $MV$ are $n$ queries, which take $O(nD^{0.437}) = O(n^{1+0.437\tau})$ time.

Conjecture 5.7 of [155]. The proof is the same, with $X := N_{[n],I}$ and $Y := V$, which are both known in Phase 2 ($X$ may have repeated columns, which do not affect the argument). Phase 1 only reads its input.

Conjecture 5.12 of [155]. Let $X := U$, an $n \times t_1$ matrix, and $Y := N_{I,J}$, a $t_1 \times t_2$ matrix, which are both known in Phase 3, with the inner dimension $D := t_1 = n^{\tau_1}$. Cut the $n$ rows of $U$ into $n/t_2$ blocks of $t_2$ rows, and preprocess each of the $n/t_2$ products of a block with $Y$, which are products of a $t_2 \times t_1$ matrix by a $t_1 \times t_2$ matrix, by Corollary 26. It applies because $t_2 = n^{\tau_2} \ge n^{18\tau_1} = D^{18}$ when $\tau_1 < \tau_2/18$. So Phase 3 costs $(n/t_2)\cdot O(t_2^2/D^{0.063}) = O(n^{1+\tau_2-0.063\tau_1})$ time. In Phase 4, $(UN_{I,J}V)_{i,j} = \sum_{\ell=1}^{t_2}(UN_{I,J})_{i,\ell}V_{\ell,j}$ is a sum of at most $t_2$ entries of $XY$, which are $t_2$ queries and take $O(t_2D^{0.437}) = O(n^{\tau_2+0.437\tau_1})$ time. Phase 2 only stores $I$.

General $\tau$. Let $0 < \tau< \varepsilon^*$, and fix $\varepsilon$ with $\tau< \varepsilon< \varepsilon^*$. By Corollary 31, with $c$ and $\theta$ chosen as in the proof of Theorem 24 for $q := \frac{1}{2}$, there is, after absorbing the logarithmic factors, a $\gamma> 0$ such that, since $D = n^\tau\le n^\varepsilon$, the arguments above go through with Phase 2 in $O(n^2/D^\gamma) = O(n^{2-\gamma\tau})$ time and Phase 3 in $O(nD^{1/2}) = O(n^{1+\tau/2})$ time. For Conjecture 5.12, we use the same blocks, which requires $n^{\tau_1} \le n^{\varepsilon\tau_2}$, that is, $\tau_1 \le\varepsilon\tau_2$ for some $\varepsilon< \varepsilon^*$.

We conclude this section by discussing which results of [155] rely on the conjectures in the refuted regime. From Conjectures 5.2 and 5.7, [155] derive trade-offs between the worst-case update time $u$ and query time $q$ of dynamic algorithms, for every $0 < \tau< 1$: either $u = \Omega(n^{\omega(1,1,\tau)-\tau-\varepsilon})$ or $q = \Omega(n^{1+\tau-\varepsilon})$. These trade-offs cover the dynamic matrix product, inverse, and adjoint problems with column updates and row queries (Theorem 5.3 and Corollary 5.5 in [155]) or with element updates and row queries (Theorem 5.8 and Corollary 5.9 in [155]), and dynamic transitive closure and DAG path counting with vertex or edge updates and source queries (Corollaries 5.4 and 5.10 in [155]). For element or pair queries, the condition becomes $q = \Omega(n^{\tau-\varepsilon})$ (Corollaries 5.9 and 5.10 in [155]).

The main bounds, such as $u+q = \Omega(n^{1.528})$ for column updates and row queries, balance the two terms at $\tau\approx0.53$. For this setting, the current best upper bound on $\omega(1,1,\tau)$ is strictly greater than 2, so the products are far from thin enough for our construction, and these bounds are still supported by the unrefuted forms of the conjectures. The same holds for the bounds from Conjecture 5.12 in [155], such as $\Omega(n^{1.407})$ for dynamic determinant, rank, and bipartite perfect matching (Corollaries 5.15 and 5.16 in [155]), since Corollary 40 refutes Conjecture 5.12 in [155] only for $\tau_1 < 0.1204\tau_2$, far from the $\tau_1+\tau_2 \approx1.41$ at which these bounds are proved.

For $\tau< 0.1204$, however, where $\omega(1,1,\tau) = 2$, the trade-offs say that an update requires $n^{2-\tau-o(1)}$ time unless a row or source query takes $n^{1+\tau-o(1)}$ time, or, in the versions with element or pair queries, unless such a query takes $n^{\tau-o(1)}$ time; Corollary 40 refutes the forms of the conjectures that would imply these trade-offs. This is the end of the trade-offs with the fastest queries, where they match the algorithms of [141], [155], such as element updates in $O(n^{2-\tau})$ time with element queries in $O(n^\tau)$ time [141], Theorem 3]. For small $\tau$, these algorithms are no longer known to be conditionally optimal.

The same goes for results that rely on this end of the trade-offs. Haeupler, Long, and Saranurak [101] cite Corollary 5.10 of [155] as showing that incremental and decremental distance oracles with approximation factor below $5/3$ cannot have worst-case update time $n^{2-\Omega(1)}$ and query time $n^{o(1)}$. With queries this fast, ruling out update time $n^{2-\delta}$ needs Conjecture 5.7 of [155] for some $\tau< \delta$, so new evidence of hardness would be needed for $\delta\le0.1204$, although the remaining cases of the conjecture still rule out update time $O(n^{1.87})$. Van den Brand, Song, and Zhou [156] show that their data structure for dynamic attention, with amortized update time $n^{\omega(1,1,\tau)-\tau}$, is conditionally optimal for every $0 < \tau\le1$ under a variant of Conjecture 5.7 of [155] in which Phase 2 gives a matrix with at most $n^\tau$ nonzero entries instead of the vector $I$. Our proof applies to this variant as well, since after Phase 2 the product again has inner dimension at most $n^\tau$, so this lower bound needs a new hypothesis for $\tau< 0.1204$.

Conclusion

We believe our results open exciting research directions in both algorithm design and complexity theory, and we highlight a few here.

Fine-grained complexity. An important message from our paper is that the fine-grained reductions that FGC has built up are incredibly valuable, and that it is even more important to design more reductions to enable further algorithm design and improved conditional hardness.

Which hypotheses should take the place of 3SUM and APSP? As discussed in the introduction, one possibility is to restrict these hypotheses to combinatorial algorithms [157, 154, 26], motivated by the search for more practical algorithms or practical hardness. Another possibility is the hypothesis that the balanced case of All-Edges Sparse Triangle, with $m$ edges, needs $m^{4/3-o(1)}$ time. Our algorithms do not speed up this balanced case, and since it is the target of several lower-bound reductions from 3SUM and Exact Triangle [136, 116, 160, 71], the hypothesis still gives many problems evidence of hardness. Is it possible that the lopsided version has a faster algorithm but the balanced version does not?

More broadly, the complexity of the many problems that are only known to be 3SUM-hard or APSP-hard, rather than equivalent to 3SUM or APSP, is now open (the gray boxes of Figure 1); we get no faster algorithms for them, since the reductions go the other way. For instance, can one decide in truly subquadratic time whether $n$ points in the plane contain three on a line [92]? As far as we know, the fastest known running time is still $O(n^2)$, without even any polylogarithmic improvements. A similar question arises for the hypotheses that our results say nothing about, such as $k$-SUM for $k \ge4$, 3SUM-Indexing, OMv (without hints), and Orthogonal Vectors. Do ideas like ours help to speed up algorithms for any of them?

Our new technique gives a way to deal with sparsity, at least in lopsided instances. Much work has gone into algorithms for sparse matrix multiplication [171, 25, 137, 2, 38, 82, 95]. However, our techniques do not seem to give new algorithms for sparse matrix multiplication. Is improving upon the known algorithms for this problem hard, or is there another technique that could help here? This question is related to the balanced version of All-Edges Sparse Triangle.

What are the true exponents of 3SUM, APSP, and the many other problems of Figure 1 that now have faster algorithms? Our exponents can certainly be improved, and one place to gain is the reductions. Each of our algorithms composes known reductions with the matrix theorem, and many reductions lose a constant fraction of the saving in the exponent. For instance, the reduction from Exact Triangle to Lopsided All-Edges Sparse Triangle keeps only half of the saving in the exponent, and the reductions from 3SUM and APSP to Exact Triangle keep a half and a third (Section 3). These reductions, like most fine-grained reductions, were designed to prove hardness, where it only matters that they keep some polynomial saving. Now that they are being used as algorithms, how much they keep is more important. Which of them can be improved? Are there more efficient reductions from all these different problems to Lopsided All-Edges Sparse Triangle, or to thin matrix products themselves, that skip the intermediate problems?¹⁰

Some of these losses cannot be removed by reductions alone. Sheffield, Vassilevska Williams, and Xi [149] show that the loss of a third in the reduction from the all-edges version of a triangle problem to its detection version [157, 159] is optimal for black-box reductions, those that work for every relation on the three edge weights at once, so a tighter reduction must use the structure of the particular problem. The loss can still be avoided on the algorithmic side, as in the footnote: our algorithm for Exact Triangle solves the all-edges version at no extra cost, as all known algorithms for triangle problems do [149]. Our reductions are not black-box in this sense either: the reduction from Exact Triangle to Lopsided All-Edges Sparse Triangle hashes the weights modulo a prime, which uses the additive structure of the problem, so the barriers of [149] say nothing against improving the half that it loses. The same paper also gives new black-box reductions among triangle, listing, and matrix problems, and such reductions are now of algorithmic interest as well.

Beyond the reductions, all of our algorithms use the same matrix theorem, and a different technique might give improvements such as better exponents, combinatorial or practical algorithms, or algorithms for All-Edges Sparse Triangle with a middle part larger than $n^{0.12}$, such as the $n^{1/2}$ that arises in an approach to girth approximation [139].

Algebraic algorithm design. Any improvement to the matrix theorem carries over to all the other speedups we achieve here. We do not know the true complexity of computing a set $W$ of entries of a thin product $XY$, as a function of $D$ and $|W|$. In particular, to handle $D > N^{0.12}$, one could look for new base identities that compute longer inner products per multiplication than Schönhage’s, perhaps by computer search [84, 114, 134], or ask whether other known identities, such as border rank expressions for the Coppersmith–Winograd tensor [69], which underlie the current bounds on $\omega$ and $\alpha$ [161, 10, 72], can be used in the same way. It is not immediately clear how to make use of other tensors, since our analysis uses specific properties of Schönhage’s identity, such as its sparsity and the way the leaves of its recursion are shared among the output entries (Section 2.4). If other tensors cannot be used in this way, it would be interesting to understand why. For instance, even if the rank or border rank of a tensor is known, it might be that we need to find different rank or border rank expressions to use, similar to work on improving leading constants of matrix multiplication algorithms (see, e.g., [36, 118, 33, 31, 131]).

Acknowledgments and Methodology

How the result was originally found and shared with the authors. An Anthropic employee used an internal research model to investigate open problems in the theory of cryptography. One of them was about cryptographic constructions based on the average-case hardness of Zero-$k$-Clique [121, 20]. Claude was tasked with verifying and improving the constructions, but instead developed this algorithm, first for the average case, then for the worst case. The session used 16M output tokens with no human input.

Anthropic shared the algorithm with the authors in September 2026 under a confidentiality agreement, offered compensation, and provided access to the public version of Claude.

Shared algorithm. The algorithm that was shared with the authors is essentially the algorithm in Section 2, although presented differently and with other numerical parameters. A different reduction from Exact Triangle to the problem of Section 2 was also shared, although the authors did not include it here, and instead used known reductions to All-Edges Sparse Triangle, which suffice and strengthen the connection. The authors derived the data structure version in Section 4 and its connection with Hinted OMv.

Lean certification. After the paper was completely written, Anthropic used an internal research model to certify this paper’s main results using the Lean 4 proof assistant with the Mathlib library. The formalization is available at https://github.com/anthropics/formal-math/tree/main/3sum-apsp. The following are formalized there: Theorem 19 (Exact Triangle in $O(n^{2.9983})$ time), Theorem 22 (3SUM in $O(n^{1.9992})$ time, and the $(\min,+)$-product and APSP in $O(n^{2.99942})$ time), and the Zero-Weight case of Corollary 39 (Zero-Weight $k$-Clique in $O(n^{k-0.0017\lfloor k/3\rfloor})$ time for every $k \ge3$). All lemmas and prior results that these statements rely on are also formalized there.

Further acknowledgments. The authors used Claude to help with writing, figure creation, and checking mathematical details throughout this project. They have attempted to make the entire paper accessible and well-written, and take full responsibility for its contents.

The authors’ contribution to this publication was done in their individual capacities, not in connection with their duties or responsibilities to Columbia University or to MIT. Any opinions expressed in this work are opinions of the authors; they are not official positions or endorsements of Columbia University or of MIT, and do not represent the research of the authors at Columbia University or MIT.