View source for Freivalds' algorithm

'''Freivalds' algorithm''' (named after [[Rūsiņš Mārtiņš Freivalds]]) is a probabilistic [[randomized algorithm]] used to verify [[matrix multiplication]]. Given three ''n''&nbsp;&times;&nbsp;''n'' [[matrix (mathematics)|matrices]] <math>A</math>, <math>B</math>, and <math>C</math>, a general problem is to verify whether <math>A \times B = C</math>. A naïve [[algorithm]] would compute the product <math>A \times B</math> explicitly and compare term by term whether this product equals <math>C</math>. However, the [[Computational complexity of matrix multiplication|best known matrix multiplication algorithm]] runs in <math>O(n^{2.3729})</math> time.<ref name="williams">{{cite CiteSeerX |title=Breaking the Coppersmith-Winograd barrier
|author=Virginia Vassilevska Williams|citeseerx=10.1.1.228.9947
|authorlink= Virginia Vassilevska Williams}}</ref> Freivalds' algorithm utilizes [[randomization]] in order to reduce this time bound to <math>O(n^2)</math><ref>
{{cite journal
|doi=10.1145/234313.234327
|last=Raghavan
|first=Prabhakar
|title=Randomized algorithms
|journal=ACM Computing Surveys
|year=1997
|volume=28
|pages=33–37
|s2cid=207196543
|accessdate=2008-12-16
|url=http://portal.acm.org/citation.cfm?id=234327|doi-access=free
}}</ref>
[[with high probability]]. In <math>O(kn^2)</math> time the algorithm can verify a matrix product with probability of failure less than <math>2^{-k}</math>.

==The algorithm==

===Input===
Three ''n''&nbsp;&times;&nbsp;''n'' [[matrix (mathematics)|matrices]] <math>A</math>, <math>B</math>, and <math>C</math>.

===Output===
Yes, if <math>A \times B = C</math>; No, otherwise.

===Procedure===
# Generate an ''n''&nbsp;&times;&nbsp;1 [[random]] 0/1 [[vector (geometric)|vector]] <math>\vec{r}</math>.
# Compute <math>\vec{P} = A \times (B \vec{r}) - C\vec{r}</math>.
# Output "Yes" if <math>\vec{P} = (0,0,\ldots,0)^T</math>; "No," otherwise.

===Error===
If <math>A \times B = C</math>, then the algorithm always returns "Yes". If <math>A \times B \neq C</math>, 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 <math>O(kn^2)</math> and error probability of <math>\leq 1/2^k</math> is achieved.

==Example==
Suppose one wished to determine whether:
:<math>
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.
</math>
A random two-element vector with entries equal to 0 or 1 is selected{{snd}} say <math>\vec{r} = \begin{bmatrix}1 \\ 1\end{bmatrix}</math>{{snd}} and used to compute:
:<math>
\begin{align}
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{align}
</math>
This yields the zero vector, suggesting the possibility that AB = C. However, if in a second trial the vector <math>\vec{r} = \begin{bmatrix}1 \\ 0\end{bmatrix}</math> is selected, the result becomes:
:<math>
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}.
</math>
The result is nonzero, proving that in fact AB &ne; C.

There are four two-element 0/1 vectors, and half of them give the zero vector in this case (<math>\vec{r} = \begin{bmatrix}0 \\ 0\end{bmatrix}</math> and <math>\vec{r} = \begin{bmatrix}1 \\ 1\end{bmatrix}</math>), so the chance of randomly selecting these in two trials (and falsely concluding that AB=C) is 1/2<sup>2</sup> 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''&nbsp;&times;&nbsp;''B'' = ''C'', then ''p'' = 0, and if ''A''&nbsp;&times;&nbsp;''B'' ≠ ''C'', then ''p'' ≤ 1/2.

===Case ''A''&nbsp;&times;&nbsp;''B'' = ''C''===
:<math>
\begin{align}
\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{align}
</math>

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

:<math>\Pr[\vec{P} \neq 0] = 0</math>

===Case ''A''&nbsp;&times;&nbsp;''B'' ≠ ''C''===
Let <math>D</math> such that

:<math>\vec{P} = D \times \vec{r} = (p_1, p_2, \dots, p_n)^T</math>

Where

:<math>D = A \times B - C = (d_{ij})</math>.

Since <math>A \times B \neq C</math>, we have that some element of <math>D</math> is nonzero. Suppose that the element <math>d_{ij} \neq 0</math>. By the definition of [[matrix multiplication]], we have:

:<math>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</math>.

For some constant <math>y</math>.
Using [[Bayes' theorem]], we can partition over <math>y</math>:

{{NumBlk|:|<math>\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]</math>|{{EquationRef|1}}}}

We use that:

:<math>\Pr[p_i = 0 | y = 0] = \Pr[r_j = 0] = \frac{1}{2}</math>
:<math>\Pr[p_i = 0 | y \neq 0] = \Pr[r_j = 1 \land d_{ij}=-y] \leq \Pr[r_j = 1] = \frac{1}{2}</math>

Plugging these in the equation ({{EquationNote|1}}), we get:

:<math>
\begin{align}
\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{align}
</math>

Therefore, 
:<math>\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}.</math>
This completes the proof.

==Ramifications==
Simple algorithmic analysis shows that the running time of this [[algorithm]] is <math>O(n^2)</math> (in [[big O notation]]). This beats the classical [[deterministic algorithm|deterministic algorithm's]] runtime of <math>O(n^3)</math> (or <math>O(n^{2.373})</math> if using [[Computational complexity of matrix multiplication|fast matrix multiplication]]). The error analysis also shows that if the [[algorithm]] is run <math>k</math> times, an [[error bound]] of less than <math>1/2^k</math> can be achieved, 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]].

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

==See also==
* [[Schwartz–Zippel lemma]]

==References==
{{reflist}}
* Freivalds, R. (1977), “Probabilistic Machines Can Use Less Running Time”, IFIP Congress 1977, pp.&nbsp;839–842.
*{{Citation
| publisher = Cambridge University Press
| last1 = Mitzenmacher
| first1 = Michael
| author1-link = Michael Mitzenmacher
| last2 = Upfal
| first2 = Eli
| author2-link = Eli Upfal
| title = Probability and computing: Randomized algorithms and probabilistic analysis
| year = 2005
| url = https://books.google.com/books?id=0bAYl6d7hvkC&pg=PA8
| pages = 8–12
| isbn = 0521835402
}}

{{Numerical linear algebra}}

{{DEFAULTSORT:Freivald's Algorithm}}
[[Category:Articles containing proofs]]
[[Category:Matrix theory]]
[[Category:Randomized algorithms]]
[[Category:Matrix multiplication algorithms]]