---
title: Point-Approximate Matrix Multiplication (PAMM)
url: https://www.emergentmind.com/topics/point-approximate-matrix-multiplication-pamm
type: topic
---

# Point-Approximate Matrix Multiplication (PAMM)

Searching arXiv for recent and foundational papers relevant to approximate and pointwise matrix multiplication.
Point-Approximate Matrix Multiplication (PAMM) denotes a family of approximate matrix multiplication formulations in which the objective is not merely to compress or sketch operands, but to approximate the product \(C=AB\) with control that is meaningful at the level of specific output points, bilinear forms, or entrywise interactions. Across the literature, however, the term is not standardized. Some papers directly support a pointwise interpretation through entrywise additive or max-norm guarantees; others are better understood as PAMM-adjacent because they approximate the full product in operator or Frobenius norm, which then induces control on entries or bilinear queries. A second nearby line of work replaces pointwise screening by hierarchical blockwise screening in the \((i,j,k)\) product space, most prominently Sparse Approximate Matrix Multiply (SpAMM), which is closely related in spirit but differs in granularity [1005.2724], [1011.3534], [1906.08148], [2504.14723].

## 1. Definitional scope and terminology

PAMM does not appear in the surveyed papers as a universally adopted formal label. The literature instead splits across several approximation notions.

One strand studies **operator-norm approximate matrix multiplication for a fixed pair of matrices**. For \(A\in\mathbb R^{n\times m}\) and \(B\in\mathbb R^{n\times p}\), the goal is to construct sketches \(\widetilde A,\widetilde B\) such that
\[
\|\widetilde A^\top \widetilde B - A^\top B\| \le \varepsilon \|A\|\,\|B\|
\]
with high probability [1005.2724]. This is highly relevant to PAMM because it is a fixed-instance guarantee and, in the random-projection case, even yields a bilinear-form statement uniform in \(x,y\):
\[
|x^\top(\widetilde A^\top \widetilde B-A^\top B)y| \le \varepsilon \|Ax\|\,\|By\|
\]
for all \(x,y\) [1005.2724].

A closely related strand studies **spectral-norm AMM via Gaussian projection** in the form
\[
XY \approx XG^\top GY,
\]
with guarantee
\[
\|XG^\top GY-XY\|_2 \le \varepsilon \|X\|_2\|Y\|_2,
\]
where the sketch dimension depends on the nuclear ranks of \(X\) and \(Y\), rather than ambient dimension [1403.7683]. This is not entrywise PAMM in the narrow sense, but it implies additive control of every bilinear query \(u^\top(XY)v\) [1403.7683].

A different strand is genuinely **entrywise or coordinatewise**. For arithmetic \(0\)-\(1\) matrix multiplication \(C=AB\), with
\[
A\in\{0,1\}^{p\times q},\qquad B\in\{0,1\}^{q\times r},
\]
the paper on Hamming-distance reductions proves that one can compute an approximation \(C''\) such that, with high probability, every entry satisfies
\[
|C''_{ij}-C_{ij}| \le \epsilon \min\{\mathrm{ham}(A_{i*},B_{*j}),\ q-\mathrm{ham}(A_{i*},B_{*j})\}
\]
[2504.14723]. This is the clearest direct PAMM result in the provided corpus.

A further line concerns **entrywise max-norm guarantees in coded distributed multiplication**. In the coded setting with \(\|A\|_F,\|B\|_F\le \eta\), the approximate recovery threshold is defined through the requirement
\[
|\widehat C_{\mathcal S}[i,j]-C[i,j]| \le \epsilon,\qquad (i,j\in[n]),
\]
for all valid worker sets \(\mathcal S\) of size at least \(K\) [2105.01973]. This is again pointwise, but in a distributed coding model.

By contrast, SpAMM and related methods are best described as **hierarchical blockwise approximate multiplication**, not pointwise in the scalar sense. They decide whether a whole submatrix product should be omitted according to a norm test in the product space [1011.3534], [1203.1692], [1906.08148], [2103.13042]. This suggests a useful distinction: PAMM in the narrow sense is scalar- or query-oriented; SpAMM is a block-hierarchical analogue.

## 2. Operator-norm and bilinear-form PAMM

The paper on low-rank matrix-valued Chernoff bounds gives an early intrinsic-dimension-dependent theory for approximate multiplication in spectral norm [1005.2724]. Its target product is \(A^\top B\), with \(A\in\mathbb R^{n\times m}\), \(B\in\mathbb R^{n\times p}\), and the objective is to produce sketches \(\widetilde A,\widetilde B\) such that
\[
\|\widetilde A^\top \widetilde B-A^\top B\| \le \varepsilon \|A\|\,\|B\|.
\]
Two mechanisms are analyzed.

The first uses a random sign matrix \(R\in\mathbb R^{t\times n}\), scaled by \(1/\sqrt t\), and sets
\[
\widetilde A=RA,\qquad \widetilde B=RB.
\]
If both matrices have rank at most \(r\), then \(t=\Omega(r/\varepsilon^2)\) suffices for the stronger bilinear guarantee
\[
|x^\top(\widetilde A^\top \widetilde B-A^\top B)y| \le \varepsilon \|Ax\|\,\|By\|
\]
for all \(x,y\), with probability at least \(1-e^{-\Omega(r)}\) [1005.2724]. If instead the stable ranks are at most \(\widetilde r\), then \(t=\Omega(\widetilde r/\varepsilon^4)\) suffices for the operator-norm guarantee with probability at least \(1-e^{-\Omega(\widetilde r/\varepsilon^2)}\) [1005.2724].

The second mechanism uses nonuniform row sampling with
\[
p_i = \frac{\|A_{(i)}\|\,\|B_{(i)}\|}{\sum_{j=1}^n \|A_{(j)}\|\,\|B_{(j)}\|}.
\]
Sampling \(t=\Omega(\widetilde r\log(\widetilde r/\varepsilon^2)/\varepsilon^2)\) i.i.d. rows and rescaling them yields
\[
\|\widetilde A^\top \widetilde B-A^\top B\| \le \varepsilon \|A\|\,\|B\|
\]
with probability \(1-1/\mathrm{poly}(\widetilde r)\) [1005.2724]. The significance is that the sketch size depends on rank or stable rank rather than the ambient dimensions \(m,p,n\).

A closely related result replaces stable rank by **nuclear rank** in the Gaussian projection model [1403.7683]. There, for \(X\in\mathbb R^{n\times d}\), \(Y\in\mathbb R^{d\times m}\), and Gaussian \(G\in\mathbb R^{t\times d}\) with \(G_{ij}\sim\mathcal N(0,1/t)\), one sets
\[
\widehat X=XG^\top,\qquad \widehat Y=GY,
\]
and obtains
\[
\|XG^\top GY-XY\|_2 \le \varepsilon \|X\|_2\|Y\|_2
\]
with probability at least \(1-\delta\), provided
\[
t=\Omega\!\left(\frac{\operatorname{nr}(X)+\operatorname{nr}(Y)+\log\log(1/\varepsilon)+\log(1/\delta)}{\varepsilon^2}\right),
\]
where
\[
\operatorname{nr}(A)=\frac{\|A\|_*}{\|A\|_2}
\]
[1403.7683]. Since a spectral-norm bound controls all bilinear forms, this gives uniform additive point-query control, even though the theorem is global rather than entrywise.

A plausible implication is that operator-norm AMM papers provide a broad PAMM foundation whenever the end task is access to bilinear forms, Gram entries, or matrix-vector queries, rather than strict per-entry relative error.

## 3. Direct entrywise PAMM for arithmetic \(0\)-\(1\) multiplication

The most explicit point-approximate formulation in the surveyed material arises from the reduction between arithmetic \(0\)-\(1\) matrix multiplication and all-pairs Hamming distances [2504.14723].

For
\[
A\in\{0,1\}^{p\times q},\qquad B\in\{0,1\}^{q\times r},
\]
the arithmetic product is
\[
C_{ij}=\sum_{\ell=1}^q A_{i\ell}B_{\ell j}.
\]
The paper proves a reverse linear-time reduction from all row-column Hamming distances to arithmetic \(0\)-\(1\) multiplication. For binary vectors \(x,y\in\{0,1\}^q\),
\[
x\cdot y=\frac{x^1+y^1-\mathrm{ham}(x,y)}{2},
\]
and therefore
\[
C_{ij}=\frac{A_{i*}^1+B_{*j}^1-\mathrm{ham}(A_{i*},B_{*j})}{2}
\]
[2504.14723]. Given all row-column Hamming distances, the exact product can thus be recovered in
\[
O(pq+qr+pr)
\]
time [2504.14723].

The approximation enters through a randomized algorithm for all-pairs Hamming distances. For \(\delta\in(0,\tfrac12)\), it computes \(\delta\)-approximations of all distances \(\mathrm{ham}(A_{i*},B_{*j})\) in time
\[
O\!\left((\log N/\delta^2)\bigl(N\log_{1+\delta}q + pr\log\log_{1+\delta}q\bigr)\right),
\qquad N=pq+qr,
\]
with high probability [2504.14723]. The algorithm is multiscale: for thresholds \(t=1,(1+\delta),(1+\delta)^2,\dots\), it constructs a randomized map
\[
f:\{0,1\}^d\to\{0,1\}^k,\qquad f(x)=Fx^\top \pmod 2,
\]
where each entry of the random matrix \(F\) is \(1\) independently with probability \(1/(4t)\) [2504.14723]. Threshold tests on projected Hamming distances distinguish whether \(\mathrm{ham}(x,y)\le t\) or \(\mathrm{ham}(x,y)\ge (1+\epsilon)t\).

To sharpen the product bound, the paper also approximates
\[
\mathrm{ham}(A_{i*},\bar B_{*j})=q-\mathrm{ham}(A_{i*},B_{*j}),
\]
where \(\bar B\) is the bitwise complement of \(B\) [2504.14723]. This yields two candidate estimators,
\[
D_{ij}=\frac{A_{i*}^1+B_{*j}^1-W_{ij}}{2},
\qquad
D'_{ij}=\frac{A_{i*}^1+B_{*j}^1-(q-W'_{ij})}{2},
\]
where \(W_{ij}\) and \(W'_{ij}\) approximate \(\mathrm{ham}(A_{i*},B_{*j})\) and \(q-\mathrm{ham}(A_{i*},B_{*j})\), respectively [2504.14723].

The final theorem states that one can compute an approximation \(C''\) such that, with high probability, for all \(i,j\),
\[
|C''_{ij}-C_{ij}| \le \epsilon \min\{\mathrm{ham}(A_{i*},B_{*j}),\ q-\mathrm{ham}(A_{i*},B_{*j})\},
\]
in time
\[
O\!\left((\log N/\epsilon^2)\bigl(N\log_{1+\epsilon}q + pr\log\log_{1+\epsilon}q\bigr)\right)
\]
[2504.14723]. This is a genuine PAMM theorem: the guarantee is per entry, not merely in a global matrix norm.

A limitation is equally clear from the statement. The error is additive and geometry-dependent; if both \(\mathrm{ham}(A_{i*},B_{*j})\) and \(q-\mathrm{ham}(A_{i*},B_{*j})\) are large, the bound can be weak. The method is also specific to arithmetic multiplication of binary matrices.

## 4. Pointwise estimators and max-norm recovery in randomized and coded settings

The 2025 mean-estimation perspective on randomized and quantum AMM is not framed as PAMM, but it contains explicit coordinatewise estimators that are directly reusable for pointwise approximation [2510.08509].

For column-sampling AMM, if \(A\in\mathbb R^{n\times q}\), \(B\in\mathbb R^{q\times m}\), and \(t\) is sampled with probability \(p_t\), then
\[
(AB)_{ij}
=
\mathbb E_{t\sim p}\!\left[\frac{A_{it}B_{tj}}{p_t}\right].
\]
The paper gives the corresponding variance formula
\[
\operatorname{Var}\!\left[\frac{A_{it}B_{tj}}{p_t}\right]
=
\sum_k \frac{A_{ik}^2B_{kj}^2}{p_k}-(AB)_{ij}^2
\]
[2510.08509]. This is a direct scalar PAMM estimator, even though the paper uses it to control the full matrix in max norm or Frobenius norm.

For the random-walk estimator for \(A_1\cdots A_k\), the paper constructs a matrix-valued random variable \(X\) with
\[
\mathbb E[X]=A_1\cdots A_k,
\]
and proves
\[
\operatorname{Var}[X_{ij}] \le \frac{(d_0)_i}{q_i}(A_1\cdots A_k)_{ij}
\]
[2510.08509]. Again, that is essentially a pointwise statement, though the theorems target full-matrix max-norm and Frobenius error.

The coded distributed setting goes further by making the **formal guarantee itself entrywise**. In \(\epsilon\)-approximate coded matrix multiplication, \(A,B\in\mathbb R^{n\times n}\) satisfy
\[
\|A\|_F\le \eta,\qquad \|B\|_F\le \eta,
\]
each worker stores a \(1/m\) fraction of each multiplicand, and the approximate recovery threshold is defined as the smallest \(K\) such that for all subsets \(\mathcal S\) with \(|\mathcal S|\ge K\),
\[
|\widehat C_{\mathcal S}[i,j]-C[i,j]|\le \epsilon
\]
for all entries \((i,j)\) [2105.01973].

The main theorem states that
\[
K^*(m,\epsilon)=m,
\]
whereas exact recovery requires \(2m-1\) workers [2105.01973]. Achievability is obtained through approximate MatDot codes: if evaluation points satisfy
\[
|\lambda_i|< \frac{\epsilon}{6\eta^2\sqrt{2m-1}(m-1)m},
\]
then recovery threshold \(m\) is achievable [2105.01973]. The converse shows \(K^*(m,\epsilon)\ge m\) for \(0<\epsilon<\eta^2\), so the threshold is optimal in that regime [2105.01973].

This is PAMM in a strong distributed sense: the approximation is full-matrix, but the metric is coordinatewise max norm. A common misconception is that approximation in coded multiplication merely relaxes average accuracy. Here the theorem is deterministic and per-entry under bounded Frobenius norms [2105.01973].

## 5. Product-space PAMM and the SpAMM lineage

SpAMM is not scalar-pointwise PAMM, but it is one of the clearest realizations of approximate multiplication by screening insignificant contributions during the multiplication itself [1011.3534], [1203.1692]. Its target setting is matrices with decay, such as
\[
|a_{ij}| < c\lambda^{|i-j|}
\]
or
\[
|a_{ij}| < \frac{c}{|i-j|^\lambda+1},
\]
where entries become small away from a locality structure [1011.3534], [1203.1692].

Matrices are stored recursively as quadtrees. At each level,
\[
A^k=\begin{pmatrix}
A_{11}^{k+1} & A_{12}^{k+1}\\
A_{21}^{k+1} & A_{22}^{k+1}
\end{pmatrix},
\]
and candidate subproducts are screened using the Frobenius norm. The defining SpAMM rule is
\[
\|A^k\|\,\|B^k\|<\tau \quad \Rightarrow \quad A^kB^k \text{ is discarded}
\]
[1011.3534]. In the optimized formulation at finer granularity,
\[
A^{t+1}_{ik} B^{t+1}_{kj} \equiv 0
\quad\text{if}\quad
\|A^{t+1}_{ik}\|_F\,\|B^{t+1}_{kj}\|_F<\tau
\]
[1203.1692].

The approximation therefore occurs in the **product space** rather than by first sparsifying \(A\) and \(B\). This matters conceptually for PAMM because it is the same basic philosophy—neglect small contributions \(A_{ik}B_{kj}\)—but at block granularity rather than scalar granularity. The 2010 SpAMM paper explicitly contrasts this with matrix-space truncation and reports that, for matched application-level error in electronic energy,
SpAMM requires fewer to far fewer floating point operations than dropping on tested quantum chemical matrices [1011.3534].

The 2019 distributed recursive framework compares three methods for matrices with exponential decay: input truncation, SpAMM-style sub-matrix-product truncation, and a hybrid of both [1906.08148]. For all three methods, it proves
\[
|[E_n]_{i,j}|=O(\tau),\qquad \forall i,j,
\]
and
\[
\|E_n\|_F =
\begin{cases}
O(n^{1/2}), & n\to\infty,\\[3pt]
O(\tau^{p/2}), & \tau\to 0,\ \forall p<2.
\end{cases}
\]
It also reports that the hybrid method reduces communication by a factor of about \(2\) on matrices from chemical systems with about \(10^6\) atoms [1906.08148].

The GPU implementation, cuSpAMM, replaces recursion by tiled kernels while preserving the block-norm admissibility test
\[
\|A_{i,k}\|_F\,\|B_{k,j}\|_F \ge \tau
\]
as the criterion for executing a tile product [2103.13042]. It introduces a valid-ratio metric,
\[
\text{valid ratio}=\sum_{i=0,j=0}^{N-1,N-1}V[i,j]/BDIM^3,
\]
measuring the fraction of candidate block products that survive pruning [2103.13042]. This suggests a systems-level analogue of PAMM granularity control: approximation decisions must be coarse enough to amortize control overhead on modern hardware.

A plausible implication is that SpAMM should be regarded as a hierarchical, blockwise PAMM relative rather than a direct pointwise method. The approximation acts on recursively defined submatrix products, but the underlying criterion is still significance screening within \((i,j,k)\) interactions.

## 6. Global AMM frameworks often cited in PAMM contexts

Several additional AMM frameworks are relevant to PAMM mainly by contrast.

The deterministic linear-system approach returns \(C'\) with
\[
\|C'-AB\|_F\le \delta
\]
in
\[
O\!\left(n^2\log\frac1\delta\right)
\]
time, for bounded-entry \(n\times n\) matrices [1408.4230]. Its guarantee is global Frobenius absolute error, not pointwise control, so it is better classified as norm-approximate matrix multiplication than PAMM.

OverSketch computes
\[
\widetilde C = ASS^TB
\]
with CountSketch-style blocks and shows that if
\[
d=\frac{2}{\epsilon\theta},
\]
then
\[
\Pr\!\left(\|AB-ASS^TB\|_F^2 \le \epsilon \|A\|_F^2\|B\|_F^2\right)\ge 1-\theta
\]
even while ignoring \(e\) stragglers among any \(z/b\) workers [1811.02653]. The guarantee is for-each Frobenius-norm AMM, not entrywise PAMM.

The coded CR-sampling literature constructs unbiased estimators of the full product from sampled block outer products. In approximate weighted \(CR\) coded multiplication, blocks are sampled with probabilities
\[
\Pi_i = \frac{\|A_i\|_F\|B_i\|_F}{\sum_{l=1}^K \|A_l\|_F\|B_l\|_F},
\]
yielding an unbiased estimator whose expected squared Frobenius error scales like
\[
O\!\left(\frac{\|A\|_F^2\|B\|_F^2}{t}\right)
\]
[2011.09709]. The distributed coded version based on MatDot reduces recovery threshold from \(2t-1\) to \(2(t/\rho)-1\) under compression by \(\rho\) [2011.09709]. Likewise, coded random sampling with MatDot recovers an unbiased approximation from any \(2s-1\) workers, with expected error
\[
\mathbb E\left[\|AB-\hat A\hat B\|_F^2\right]
=
\frac{1}{s}\left(\sum_{q=0}^{m-1}\|A_qB_q\|_F\right)^2
-\frac{1}{s}\|AB\|_F^2
\]
for independent sampling [1905.06942]. These are globally normed AMM estimators, though their summand-wise structure is often PAMM-adjacent.

Finally, decomposition-based approximate multiplication via TSVD, circulant decomposition, Fourier sparsification, and cycle decomposition uses the first-order expansion
\[
AB \approx A_{\text{keep}}B + \Delta A\,B_{\text{keep}},
\qquad
AB-\tilde C = \Delta A\,\Delta B,
\]
with reported \(\mathcal O(n^2\log n)\) arithmetic complexity for \(n\times n\) matrices at usable relative error around \(1\%\) [2504.19308]. This is again global Frobenius-style AMM, not pointwise PAMM, but the exact residual identity is structurally useful.

## 7. Conceptual distinctions, misconceptions, and open directions

A recurring misconception is that “approximate matrix multiplication” is a single notion. The literature instead separates at least four distinct regimes.

| Regime | Typical guarantee | Representative papers |
|---|---|---|
| Entrywise / max-norm | \(|\widehat C_{ij}-C_{ij}|\le \epsilon\) or geometry-dependent per-entry error | [2504.14723], [2105.01973] |
| Bilinear / operator norm | \(\|\widetilde A^\top \widetilde B-A^\top B\|\le \varepsilon \|A\|\|B\|\) | [1005.2724], [1403.7683] |
| Frobenius AMM | \(\|AB-\tilde C\|_F\) or expected squared Frobenius error | [1408.4230], [2011.09709], [1811.02653], [1905.06942] |
| Product-space hierarchical | norm-screened omission of submatrix products | [1011.3534], [1203.1692], [1906.08148], [2103.13042] |

A second misconception is that blockwise methods such as SpAMM are merely input-sparsification schemes. They are not: the defining approximation mechanism is omission of subproducts based on
\[
\|A_{\alpha\beta}\|_F\,\|B_{\beta\gamma}\|_F<\tau,
\]
which acts directly in the multiplication process [1906.08148].

A third misconception is that approximate coded multiplication only yields average or probabilistic output error. In fact, the approximate MatDot framework proves deterministic coordinatewise recovery guarantees under norm-bounded inputs [2105.01973].

The surveyed corpus also reveals several open directions. The 2010 and 2019 SpAMM papers explicitly treat rigorous worst-case global bounds and asymptotic behavior for slow decay as unresolved or only partially resolved [1011.3534], [1906.08148]. The binary Hamming-distance reduction suggests a precise path to PAMM for discrete domains, but this does not automatically extend to general real matrices [2504.14723]. The mean-estimation perspective suggests that many full-matrix AMM estimators can be specialized entrywise, but those specializations are usually not developed as standalone point-query algorithms [2510.08509]. This suggests that a unifying modern PAMM theory would likely need to formalize the relationship between coordinatewise estimators, global norm guarantees, and hierarchical product-space screening.

In that sense, PAMM is best understood not as a single established algorithm, but as a technically coherent research direction spanning entrywise additive approximation, max-norm recovery, fixed-instance bilinear control, and hierarchical significance screening inside matrix multiplication itself [1005.2724], [2504.14723].

Source: https://www.emergentmind.com/topics/point-approximate-matrix-multiplication-pamm