Freivalds' algorithm

Freivalds' algorithm (named after Rūsiņš Mārtiņš Freivalds) is a probabilistic randomized algorithm used to verify matrix multiplication. Given three n × n matrices ${\displaystyle A}$, ${\displaystyle B}$, and ${\displaystyle C}$, a general problem is to verify whether ${\displaystyle A\times B=C}$. A naïve algorithm would compute the product ${\displaystyle A\times B}$ explicitly and compare term by term whether this product equals ${\displaystyle C}$. However, the best known matrix multiplication algorithm runs in ${\displaystyle O(n^{2.3729})}$ time.[1] Freivalds' algorithm utilizes randomization in order to reduce this time bound to ${\displaystyle O(n^{2})}$[2] with high probability. In ${\displaystyle O(kn^{2})}$ time the algorithm can verify a matrix product with probability of failure less than ${\displaystyle 2^{-k}}$.

The algorithm

Input

Three n × n matrices ${\displaystyle A}$, ${\displaystyle B}$, and ${\displaystyle C}$.

Output

Yes, if ${\displaystyle A\times B=C}$; No, otherwise.

Procedure

1. Generate an n × 1 random 0/1 vector ${\displaystyle {\vec {r}}}$.
2. Compute ${\displaystyle {\vec {P}}=A\times (B{\vec {r}})-C{\vec {r}}}$.
3. Output "Yes" if ${\displaystyle {\vec {P}}=(0,0,\ldots ,0)^{T}}$; "No," otherwise.

Error

If ${\displaystyle A\times B=C}$, then the algorithm always returns "Yes". If ${\displaystyle A\times B\neq C}$, then the probability that the algorithm returns "Yes" is less than or equal to one half. This is called one-sided error.

By iterating the algorithm k times and returning "Yes" only if all iterations yield "Yes", a runtime of ${\displaystyle O(kn^{2})}$ and error probability of ${\displaystyle \leq 1/2^{k}}$ is achieved.

Example

Suppose one wished to determine whether:

${\displaystyle AB={\begin{bmatrix}2&3\\3&4\end{bmatrix}}{\begin{bmatrix}1&0\\1&2\end{bmatrix}}{\stackrel {?}{=}}{\begin{bmatrix}6&5\\8&7\end{bmatrix}}=C.}$

A random two-element vector with entries equal to 0 or 1 is selected – say ${\displaystyle {\vec {r}}={\begin{bmatrix}1\\1\end{bmatrix}}}$ – and used to compute:

${\displaystyle {\begin{aligned}A\times (B{\vec {r}})-C{\vec {r}}&={\begin{bmatrix}2&3\\3&4\end{bmatrix}}\left({\begin{bmatrix}1&0\\1&2\end{bmatrix}}{\begin{bmatrix}1\\1\end{bmatrix}}\right)-{\begin{bmatrix}6&5\\8&7\end{bmatrix}}{\begin{bmatrix}1\\1\end{bmatrix}}\\&={\begin{bmatrix}2&3\\3&4\end{bmatrix}}{\begin{bmatrix}1\\3\end{bmatrix}}-{\begin{bmatrix}11\\15\end{bmatrix}}\\&={\begin{bmatrix}11\\15\end{bmatrix}}-{\begin{bmatrix}11\\15\end{bmatrix}}\\&={\begin{bmatrix}0\\0\end{bmatrix}}.\end{aligned}}}$

This yields the zero vector, suggesting the possibility that AB = C. However, if in a second trial the vector ${\displaystyle {\vec {r}}={\begin{bmatrix}1\\0\end{bmatrix}}}$ is selected, the result becomes:

${\displaystyle A\times (B{\vec {r}})-C{\vec {r}}={\begin{bmatrix}2&3\\3&4\end{bmatrix}}\left({\begin{bmatrix}1&0\\1&2\end{bmatrix}}{\begin{bmatrix}1\\0\end{bmatrix}}\right)-{\begin{bmatrix}6&5\\8&7\end{bmatrix}}{\begin{bmatrix}1\\0\end{bmatrix}}={\begin{bmatrix}-1\\-1\end{bmatrix}}.}$

The result is nonzero, proving that in fact AB ≠ C.

There are four two-element 0/1 vectors, and half of them give the zero vector in this case (${\displaystyle {\vec {r}}={\begin{bmatrix}0\\0\end{bmatrix}}}$ and ${\displaystyle {\vec {r}}={\begin{bmatrix}1\\1\end{bmatrix}}}$), so the chance of randomly selecting these in two trials (and falsely concluding that AB=C) is 1/2² or 1/4. In the general case, the proportion of r yielding the zero vector may be less than 1/2, and a larger number of trials (such as 20) would be used, rendering the probability of error very small.

Error analysis

Let p equal the probability of error. We claim that if A × B = C, then p = 0, and if A × B ≠ C, then p ≤ 1/2.

Case A × B = C

${\displaystyle {\begin{aligned}{\vec {P}}&=A\times (B{\vec {r}})-C{\vec {r}}\\&=(A\times B){\vec {r}}-C{\vec {r}}\\&=(A\times B-C){\vec {r}}\\&={\vec {0}}\end{aligned}}}$

This is regardless of the value of ${\displaystyle {\vec {r}}}$, since it uses only that ${\displaystyle A\times B-C=0}$. Hence the probability for error in this case is:

${\displaystyle \Pr[{\vec {P}}\neq 0]=0}$

Case A × B ≠ C

Let ${\displaystyle D}$ such that

${\displaystyle {\vec {P}}=D\times {\vec {r}}=(p_{1},p_{2},\dots ,p_{n})^{T}}$

Where

${\displaystyle D=A\times B-C=(d_{ij})}$.

Since ${\displaystyle A\times B\neq C}$, we have that some element of ${\displaystyle D}$ is nonzero. Suppose that the element ${\displaystyle d_{ij}\neq 0}$. By the definition of matrix multiplication, we have:

${\displaystyle p_{i}=\sum _{k=1}^{n}d_{ik}r_{k}=d_{i1}r_{1}+\cdots +d_{ij}r_{j}+\cdots +d_{in}r_{n}=d_{ij}r_{j}+y}$.

For some constant ${\displaystyle y}$. Using Bayes' theorem, we can partition over ${\displaystyle y}$:

${\displaystyle \Pr[p_{i}=0]=\Pr[p_{i}=0|y=0]\cdot \Pr[y=0]\,+\,\Pr[p_{i}=0|y\neq 0]\cdot \Pr[y\neq 0]}$

(1)

We use that:

${\displaystyle \Pr[p_{i}=0|y=0]=\Pr[r_{j}=0]={\frac {1}{2}}}$

${\displaystyle \Pr[p_{i}=0|y\neq 0]=\Pr[r_{j}=1\land d_{ij}=-y]\leq \Pr[r_{j}=1]={\frac {1}{2}}}$

Plugging these in the equation (1), we get:

${\displaystyle {\begin{aligned}\Pr[p_{i}=0]&\leq {\frac {1}{2}}\cdot \Pr[y=0]+{\frac {1}{2}}\cdot \Pr[y\neq 0]\\&={\frac {1}{2}}\cdot \Pr[y=0]+{\frac {1}{2}}\cdot (1-\Pr[y=0])\\&={\frac {1}{2}}\end{aligned}}}$

Therefore,

${\displaystyle \Pr[{\vec {P}}=0]=\Pr[p_{1}=0\land \dots \land p_{i}=0\land \dots \land p_{n}=0]\leq \Pr[p_{i}=0]\leq {\frac {1}{2}}.}$

This completes the proof.

Ramifications

Simple algorithmic analysis shows that the running time of this algorithm is O(n²), beating the classical deterministic algorithm's bound of O(n³). The error analysis also shows that if we run our algorithm k times, we can achieve an error bound of less than ${\displaystyle {\frac {1}{2^{k}}}}$, an exponentially small quantity. The algorithm is also fast in practice due to wide availability of fast implementations for matrix-vector products. Therefore, utilization of randomized algorithms can speed up a very slow deterministic algorithm. In fact, the best known deterministic matrix multiplication algorithm known at the current time is a variant of the Coppersmith–Winograd algorithm with an asymptotic running time of O(n^2.3729).[1]

Freivalds' algorithm frequently arises in introductions to probabilistic algorithms because of its simplicity and how it illustrates the superiority of probabilistic algorithms in practice for some problems.

See also

• Schwartz–Zippel lemma

References

1. Virginia Vassilevska Williams. "Breaking the Coppersmith-Winograd barrier" (PDF).
2. Raghavan, Prabhakar (1997). "Randomized algorithms". ACM Computing Surveys. 28: 33. doi:10.1145/234313.234327. Retrieved 2008-12-16.

• Freivalds, R. (1977), “Probabilistic Machines Can Use Less Running Time”, IFIP Congress 1977, pp. 839–842.