Lesson 4

Randomized algorithm for matrix multiplication verification

In the matrix multiplication problem, the input consists of two matrices \(A, B \in \mathbb{Z}^{n \times n}\), and the goal is to compute their product, that is, the matrix \(C \in \mathbb{Z}^{n \times n}\) such that \[C[i,j]=\sum_{k=1}^{n}A[i,k]\cdot B[k,j] \ .\] Matrix multiplication is widely used in numerical methods and also has numerous applications in other areas (for example, pattern recognition and graph algorithms). For this reason, algorithms for this problem are studied intensively.

A straightforward algorithm for multiplying matrices has running time \(O(n^{3})\): one needs to compute \(n^{2}\) entries, and each entry takes time \(O(n)\) to compute. The first improvement of this algorithm was given by Volker Strassen in 1969. His paper showed how to solve systems of linear equations faster than in cubic time. At the time, this result surprised many: it seemed that no algorithm faster than cubic time could exist. Strassen's algorithm is based on the “divide and conquer” method and has running time \(O(n^{2.81})\). The algorithm works as follows: each of the two input matrices is divided into four submatrices of size \(n/2 \times n/2\); seven products of linear combinations of these matrices are computed; the final result is assembled from these products; since seven recursive multiplications are needed, we obtain the recurrence \(T(n) \le 7T(n/2)+O(n^{2})\), which implies the upper bound \(T(n)=O(n^{\log_2 7})\). The constant \(\log_{2} 7\) has since been improved multiple times. In general, due to the importance of this problem in various areas, the following notation is standard: the exponent of matrix multiplication is the number \(\omega\) such that for every \(\varepsilon>0\) there exists an algorithm multiplying \(n \times n\) matrices in time \(O(n^{\omega + \varepsilon})\). Strassen thus showed that \(\omega \le \log_{2} 7\). The latest improvement was obtained by Alman, Duan, Vassilevska Williams, Xu, Xu and Zhou in 2024: they proved that \(\omega \le 2{,}371339\), improving the fourth (!) digit after the decimal point.

This result is purely theoretical progress: the hidden constants in such algorithms are so enormous that the algorithms themselves are completely impractical (they can outperform the naive cubic algorithm only for such huge \(n\) that \(n \times n\) matrices would not even fit into memory). At present, many believe that \(\omega=2\)—that matrices can be multiplied in time \(O(n^{2+o(1)})\), although we still do not know how to do this. There are already research papers that analyze the running times of algorithms for various problems under the assumption that \(\omega=2\).

Suppose we have at our disposal a fast matrix multiplication algorithm and we want to verify whether it correctly multiplies huge matrices. Since the matrices are enormous, we do not have enough time to compute their product ourselves using the naive cubic algorithm. It turns out that there is a simple randomized algorithm with running time \(O(n^{2})\) that verifies the equality with high probability.

Below we present Freivalds' algorithm for this task. It is a fingerprinting algorithm: instead of comparing the objects themselves (which may be expensive or difficult), we compare the value of some function on them. We expect that this function is relatively easy to compute. We then use a simple observation: if \(f(x) \neq f(y)\), then \(x \neq y\).

For the problem of verifying the result of matrix multiplication, such a function will be multiplication by a vector. Namely, take a random vector \(r\in\{0,1\}^{n}\) and check whether \[(AB)r=Cr \ .\] It turns out that in order to compute \((AB)r\), we do not need to compute \(AB\) at all! Indeed, we can first compute \(Br\), and then multiply \(A\) by \(Br\): \[(AB)r=A(Br) \ .\] Thus, both \((AB)r\) and \(Cr\) can be computed in time \(O(n^{2})\).

Theorem. Let \(A,B,C \in \mathbb{F}^{n \times n}\) be square matrices over a field \(\mathbb{F}\), and let \(r \in \{0,1\}^{n}\) be a random vector. If \(AB \neq C\), then \[\Pr[(AB)r \neq Cr] \ge \frac{1}{2}\ .\]

Proof. Let \(D=AB-C\). Since \(D\) is a nonzero matrix, it contains a nonzero entry. Without loss of generality, suppose there is a nonzero entry in the first row \(d=D[1]\), and moreover it is the first entry in that row: \(d_{1} \neq 0\). We will show that even the probability that \(dr \neq 0\) is at least \(1/2\). Suppose that all coordinates of the vector \(r\), except the first, are already fixed. Then \[dr=0 \Leftrightarrow r_{1} = -\frac{\sum_{i=2}^{n}d_ir_i}{d_1}\ .\] Thus, for fixed \(r_{2}, \dotsc, r_{n}\), at most one value of \(r_{1}\) works, and therefore the fraction of all bad \(r\) is at most one half.◼

Problem. Is this argument correct?

In Freivalds' algorithm, one can use any two other different constants instead of zero and one: the analysis stays the same.