Faster Walsh-Hadamard Transform and Matrix Multiplication over Finite Fields using Lookup Tables

When $N$ is a power of $2$, the Walsh-Hadamard transform is defined recursively by $H_{2}=\begin{bmatrix}1&1\\ 1&-1\\ \end{bmatrix},$ and

A common task in many areas of computation is to compute the length-$N$ Walsh-Hadamard transform, i.e., given as input a length-$N$ vector $v\in\mathbb{F}^{N}$, compute the vector $H_{N}v$. The most straightforward algorithm would compute this using $O(N^{2})$ arithmetic operations, but the fast Walsh-Hadamard transform (FWHT) algorithm can compute this using only $O(N\log N)$ arithmetic operations. It is widely believed that $\Omega(N\log N)$ arithmetic operations are necessary, and a substantial amount of work has gone into proving this in restricted arithmetic models²²2For instance, this is known to hold for arithmetic circuits with ‘bounded coefficients’; see e.g., [Lok09, Section 3.3]., and studying conjectures which would imply this³³3For instance, the now-refuted conjecture that the Walsh-Hadamard transform is rigid [AW17].; see, for instance, the survey [Lok09].

A natural question arises: why restrict ourselves to arithmetic models of computation? The Walsh-Hadamard transform is commonly used in practice, and speedups using non-arithmetic operations could be impactful. Nonetheless, there has not been much work on non-arithmetic algorithms for the Walsh-Hadamard transform.

A related problem is matrix multiplication: given as input two $n\times n$ matrices, compute their product. The asymptotically fastest known algorithm is algebraic, and uses $O(n^{2.373})$ arithmetic operations. That said, non-algebraic algorithmic techniques, especially lookup tables, have been used to design more practical, ‘combinatorial’ algorithms for a variant on this problem called Boolean matrix multiplication⁴⁴4In Boolean matrix multiplication, given two matrices whose entries are from $\{0,1\}$, our goal is to multiply them over the AND-OR semiring, or equivalently, to determine which entries of their product over $\mathbb{R}$ are nonzero (but not necessarily what nonzero values they take on). since the work of [ADKF70]. These techniques save logarithmic factors over the straightforward $O(n^{3})$ time algorithm – the best algorithm along these lines runs in $n^{3}\cdot\frac{(\log\log n)^{O(1)}}{(\log n)^{3}}$ time [Yu18] in the standard word RAM model – but they are considered more practical than the algebraic algorithms (which save polynomial factors in $n$). Lookup tables have also been used in practice to approximately multiply matrices [JPK⁺20, BG21].

In this paper, we show that lookup tables can be used to speed up the asymptotically fastest known algorithms for both the Walsh-Hadamard transform and (exact) matrix multiplication over finite fields. We will show, for instance, that only $o(N\log N)$ bit operations suffice to compute the length-$N$ Walsh-Hadamard transform over finite fields when we augment the arithmetic model to allow for lookup tables. This may help to explain the difficulty of proving lower bounds for this problem, and help to guide future work on arithmetic circuit lower bounds (since any lower bounds technique would need to fail when lookup tables are allowed).

We focus here on constant-sized finite fields, though our algorithms generalize to larger finite fields as well. Our algorithms are simple modifications of the usual recursive algorithms for solving these problems, and we describe below how they compare favorably to algorithms which are used in practice.

Let $\mathbb{F}_{q}$ be the finite field of size $q$, and let $N$ be a power of $2$. In this section, we give an algorithm for computing the length-$N$ Walsh-Hadamard transform, $H_{N}$, over $\mathbb{F}_{q}$. The key idea behind our new algorithm is to pick a $K=O(\log_{q}N)$ such that $q^{K}\ll N$, and first create a lookup table of the length-$K$ Walsh-Hadamard transforms of all vectors $v\in\mathbb{F}_{q}^{K}$. We will then use this lookup table in conjunction with the following standard recursive approach for computing Kronecker powers (sometimes called Yates’ algorithm), which we will apply with $M=H_{K}$:

We can now give our main algorithm for the Walsh-Hadamard transform:

Fast algebraic algorithms for matrix multiplication over a field $\mathbb{F}$ critically rely on algebraic identities which take the following form, for positive integers $q,r$, formal variables $X_{i,j},Y_{j,k},Z_{i,k}$ for $i,j,k\in[t]$, and field coefficients $\alpha_{i,j,\ell},\beta_{j,k,\ell},\gamma_{i,k,\ell}\in\mathbb{F}$ for $i,j,k\in[t]$ and $\ell\in[r]$:

As we will see next, an identity (1) can be used to design a matrix multiplication algorithm which runs using only $O(n^{\log_{t}(r)})$ field operations. For instance, Strassen’s algorithm [Str69] gives an identity (1) with $t=2$ and $r=7$, yielding exponent $\log_{2}(7)<2.81$. The matrix multiplication exponent $\omega$ is defined as the infimum of all numbers such that, for every $\varepsilon>0$, there is a sufficiently large $t$ for which one gets an identity (1) with $r\leq t^{\omega+\varepsilon}$. Indeed, it is known [Str73] that any algebraic algorithm for matrix multiplication can be converted into an identity (1) yielding the same running time in this way, up to a constant factor, so $\omega$ captures the best exponent from any possible algebraic algorithm. The fastest known matrix multiplication algorithm [CW90, AW21] shows that $\omega<2.37286$.

The standard recursive algorithm for multiplying matrices using identity (1) works as described in Algorithm 1 below. It recurses until it gets to a base case of multiplying $n\times n$ matrices when $n\leq S$ for some parameter $S$. In the usual algorithm, one picks $S$ to be a constant, so that such matrices can be multiplied in a constant number of operations; in our improvement, we will pick a larger $S$ and multiply such matrices using a lookup table.

Normally one would analyze the recursion for the number of operations performed by this algorithm when $S$ is a constant, and conclude a total operation count of $O(n^{\log_{t}(r)})$. Here we will improve on this over finite fields using a lookup table:

We showed that for two important open problems – determining the complexity of the Walsh-Hadamard transform, and determining the leading constant of Strassen’s algorithm for matrix multiplication – asymptotic improvements over the conjectured optimal arithmetic bounds are possible if one is allowed to use bit operations rather than just arithmetic operations.

Our algorithms only made use of arithmetic operations and lookup table operations, so they could extend to other models of computation as well. One natural question is whether they extend to the standard word RAM model with word size $w=O(\log n)$ for input size $n$. Indeed, operations for the lookup tables we use, (with keys and values of $O(\log n)$ bits) require only $O(1)$ word operations in this model.

It is not hard to see that our algorithm for matrix multiplication can be implemented to take advantage of this model, improving Theorem 4.2 to running time $O(n^{\tau}/(\log n)^{\tau/2})$ (improving by a factor of $\log n$).

On the other hand, our algorithm for the Walsh-Hadamard transform seemingly cannot be implemented in the word RAM model to get asymptotic savings. The main culprit is step 2 in the proof of Lemma 3.1: if we want to efficiently use the lookup table to compute $y_{j}=Mx_{j}$, then we first have to permute the bits of the $u_{i}$ vectors so that each $x_{j}$ fits in a single word. In other words, we are given the $u_{i}$ vectors in order, and would like to permute their entries to get the $x_{j}$ vectors in order instead.

Let $s=\Theta(d^{m}/w)$ be the number of words that our input fits into. In general it is known that performing a fixed permutation of the bits of a string contained in $s$ words, for $s\geq w$, requires $\Theta(s\log s)$ word operations [BMM97]. However, our particular permutation can be broken up into $s/w$ different permutations on $w$ words (e.g., the first words in the descriptions of $x_{1},\ldots,x_{w}$ are a permutation of the first words in the descriptions of $u_{1},\ldots,u_{w}$). It can thus be performed in only $O(s\log w)$ word operations.

Since this must be done at all $m$ levels of recursion, it incurs an additional $\Theta(m\cdot d^{m}\cdot\frac{\log w}{w})$ word operations. With the parameter setting we ultimately use in Theorem 3.2, this is $\Theta(N)$ word operations to perform all these permutations. By comparison, it is not too difficult to implement the usual fast Walsh-Hadamard transform to use a total of $O(N)$ word operations as well. Hence, the time to perform permutations (which doesn’t come into play in the arithmetic or bit operation models) swamps any other computational savings in this model, and another approach is needed for a speedup.

I would like to thank Dylan McKay and Ryan Williams for invaluable discussions throughout this project, and anonymous reviewers for helpful comments. This research was supported in part by a grant from the Simons Foundation (Grant Number 825870 JA).