Mathematicians still don't know the fastest way to multiply numbers
scientificamerican.com227 points by beardyw 7 days ago
227 points by beardyw 7 days ago
Back in 2024, I was trying to optimize PostgreSQL's NUMERIC data type, which is base-10000, using Karatsuba. The problem of finding the optimal threshold of when to switch to Karatsuba turned out to be really hard, since it depends on the size of both factors combined. After some hundreds of hours, I gave up, and started thinking about if there could be a simpler solution. I came to think about another idea I'd had before but abandoned, about 64-bit modernizing the digit base from 10k to 100M, but that would be a challenge due to existing data on disk. Desperate of finding a solution, I wondered if it could be fast enough to do on-the-fly conversion back and forth between base-10k and base-100M, and then realized that, yes, of course, it will be fast already for quite small N (testing shows already between 3-6 base digits). The trick basically reduced the N in O(N^2) into half, i.e. O((N/2)^2), with some O(2*N) cost for the conversion back and forth.
I had a lot of fun hacking on this idea together with the maintainer of the NUMERIC data type, and after two months the patch finally was ready and got committed:
https://git.postgresql.org/gitweb/?p=postgresql.git;a=commit...
A bit tangential, but the folks behind the GNU Multiple Precision Library (GMPLib) have the problem of choosing algorithms more or less fleshed out. They've got some fairly approachable manual pages[1] for the various algorithms they use as operand sizes scale up, where Karatsuba is only the second of six options in terms of operational complexity.
This really demonstrates the utility of pulling in a library written by experts in the field the library handles.
This is a problem I have a lot in modern programming ecosystems: How do I tell the difference between a library written by a team of experts who have spent decades optimizing everything to do with the task and a library written by one guy that's an unnecessary straightforward wrapper over the obvious implementation?
Or increasingly, a library written by an LLM referencing a pile of such “one guy” projects of varying levels of suck (from actually good to good-got-it-is-full-of-suck).
What I have been seeing recently is having a great LLM rewrite the library leads to fewer bugs and also optimized for your use case. More expensive for sure than pulling something like epub.js but when the library has infinite open issues it can be a lot better.
You defer to the advice of experts you trust. Which somehow have become harder to come by in terms of signal to noise than 20-30 years ago.
Here is the full pgsql-hackers mailing list thread where you can follow our work from initial idea to commit: https://www.postgresql.org/message-id/flat/9d8a4a42-c354-41f...
It's complicated. :-)
There is a nice picture of the "best" choice for different ranges of sizes of numbers to be multiplied at http://gmplib.org/devel/log.i7.1024.png
More context and explanation can be found at: http://gmplib.org/devel/
I have actually had a ton of success using Strassen matrix multiplication kernels with extra structure in custom CUDA kernels (e.g. a covariance matrix is symmetric positive definite, or can be represented with Cholesky, and that comes up in a ton of useful computation). It's been a couple of years, but IIRC I would find it would start to win over the standard kernels at ~n>2500 or something (and in addition to Strassen was also exploiting the explicit structural constraints of the matrix, so not a completely fair comparison).
If you’re interested, I found https://arxiv.org/abs/2505.09814v1 to beat Strassen for medium-sized and larger covariance matrices. YMMV of course. Takes a little adjustment for XX^H but it’s not so bad.
For anyone interested in digging, check out Karatsuba's algorithm (https://en.wikipedia.org/wiki/Karatsuba_algorithm), Strassen matrix multiplication (https://en.wikipedia.org/wiki/Strassen_algorithm), and Toom-Cook multiplication (https://en.wikipedia.org/wiki/Toom–Cook_multiplication).
Of course, there are also implementation considerations. For example you can speed up Strassen by recursively breaking down the matrix into sub-matrices in parallel, but only down to a point - once the sub-matrices get small enough, it becomes faster to simply do a straight Strassen computation. And it depends on your hardware. For something seemingly so simple, you can go pretty far down a rabbit hole!
> once the sub-matrices get small enough,
the next rabbit hole starts when you try to start lifting F_2 into Z for a given <n,n,n> tensor.
Why do we make computers multiply single digit numbers, instead of taking the result from a lookup table, like humans do? To answer my own question, I am assuming it would be because multiplying would still be faster than reading from a lookup table? Any ideas?
In a sense, they do exactly that! But since there are only two single-digit numbers in binary, it makes for a pretty short table.
Hardware multipliers often use a sort of base-4-ish lookup table trick as well, using the Booth-Wallace algorithm. Booth's idea is to rewrite one of the inputs in base (usually) "4", except that the digits go from -2 to +2 instead of 0 to 3. (That's five possible digits! This helps the rewriting stage not have to propagate carries. Carry propagation is very expensive.) You can use Booth in a base higher than 4, especially if you know one of the multiplicands before the other, but you run into tradeoffs pretty quickly.
Then for each digit, you select between the other input multiplied by 0 (all zeros), +1 (identity), +2 (shift left by one bit), or -1 or -2 (flip all the bits of +1 or +2, plus a correction). Since a number has about half as many digits in base 4 as in base 2, you have about half as many digits to sum as if you'd done this in base 2.
Then you sum up all those results, but since carry propagation is expensive, you mostly use "compressors", e.g. you sum up three intermediates at a time, but you do it bit-by-bit, where three 1-bit numbers add up to a 2-bit number (from 0 to 3). This is called a Wallace Tree. The point is that you are generating carries, but you aren't propagating them, just adding them back into the set of things to be summed.
At the end of the tree step, you have just two numbers left, and you add them conventionally. That's the only step that needs full carry propagation.
If you are implementing a multiply-add, or multiplying several numbers and adding up all the results or similar, then you usually only need one full carry propagation stage.
The overall circuit has quadratic area but only a logarithmic depth in gates. IIRC whether to do Booth or not is a tradeoff: at least in some circumstances the rewrite steps make it slower but smaller. Hardware tool vendors have done a lot of work to tune these circuits very tightly, using e.g. specialized gates like AOI, heuristics for how to set up the tree, etc.
If carry propagation is so expensive, why are the mathematicians, like Karatsuba, ignoring it? It seems like we need a better complexity measure.
> Booth's idea is to rewrite one of the inputs in base (usually) "4", except that the digits go from -2 to +2 instead of 0 to 3.
That's base 5 then. It needs to go from -2 to +1 if you want base 4
It's still a modified base 4, because the significance of the i'th digit is 4^i, not 5^i.
Edited to add: I'm also not sure whether real-life implementations have -0 as an option. Of course -0 could be normalized to +0, but it might be cheaper not to bother if the sign is applied after the digit selection.
No, balanced ternary, for example, uses {-1, 0, 1}. The system you're discussing is balanced quinary (base 5).
It isn't balanced quinary, but rather redundant balanced quaternary (base 4). In balanced quinary (base 5), each digit has 5x the significance of the previous one, but in Booth's encoding algorithm it's 4x.
If digit i has significance b^i, then b (the base of the exponentiation) is the base (or radix) of the number system.
The page you linked explicitly mentions the binary version of Booth encoding as having base b=2 and three signed digits {-1, 0, 1}. The quaternary version similarly has b=4 and five signed digits {-2, -1, 0, 1, 2} ... and possibly sometimes -0 in practice, not sure.
Every combinational logic function can be implemented as a look up table, which can be implemented as a hardware ROM. It is common to consider speed-power-area tradeoffs to find good implementations.
> I am assuming it would be because multiplying would still be faster than reading from a lookup table?
Even if it isn’t, it still would be a lot cheaper. With 32-bit integers, the lookup table would have 2⁶⁴ 64-bit values. If my math six right, that is 128 exabytes of read-only memory. With 64-but integers, it truly would be impractical, at 2¹²⁸ 128-bit values.
Printing it on paper and manually looking it up is probably cheaper with today’s ram-prices
Assuming the ink is free, you use some sort of scheme to take advantage of the full printable ascii character set, and you don't bother presenting key values (if necessary, you could print a key offset at the top-left of each page), just the table outputs, I'm getting a break-even point at each character being around 0.2 square nanometers.
Does the article just end after describing the problem for me only? I am left wanting for more.
if you want a little more in depth explanation of the whole multiplication thing Ican recommend TAOCP volume 2. it has a section called 'How fast can we multiply?' it should provide more insight. there is a paper from djb (Daniel J Bernstein),which I can also recommemd: https://cr.yp.to/lineartime/multapps-20080515.pdf
Amazed I hadn't heard of this before. Would be interesting to see if they can prove that they have discovered the fastest at O(n × log n) or whether there is more still to come.
I was under the impression that because the grade school technique we learn is really just convolution over the digits, the fastest algorithms achieve o(n logn) via fourier transforms. Is that not the case?
Can you explain some more? Thanks.
If you try multiplying two numbers, say 4 digits each, using the standard grade school technique, you'll see that you end up writing the digit-wise product of each number, with a shift (adding a zero at the beginning) for each digit, followed by a sum. This is literally just a convolution operation over the digits (this is because secretly writing a number in base B is equivalent to expressing a polynomial evaluated at the integer B, and multiplying polynomials is a convolution over the coefficients). By the convolution theorem, this O(N^2) time operation can therefore be accomplished in O(N logN) time by doing multiplication in frequency space and then transforming back. This is because the Fourier transform is O(N logN), via the FFT, and multiplication is O(N).
I just looked up the answer to my original question - the Fourier trick is notionally only O(N logN), but because the FFT takes you from integers to floating points, as N gets larger you need to encode more and more bits to achieve enough precision to yield absolute errors <1 after doing both Fourier transforms. The need to encode those extra bits tacks on another logN, taking you to O(N log^2 N).
Jesus, I took 4 years of Electrical Engineering and at no point did any professor make this brilliant analogy. It would have helped many students, including myself.
I don't know what age range "grade school" is, but I remember being taught that method when I was about 7 or 8, although it didn't really "land" properly until I read the short story "The Feeling of Power" by Isaac Asimov.
What I'm surprised to see left out here (unless I missed it in the page's horrible formatting) is a mention of the way that computers multiply two integers. They use a technique I saw described in a book when I was about 11 as the "Russian Farmer Method" (or something like that, it was in English and I might have misremembered it).
In that you shift the multiplier right and multiplicand left, halving one and doubling the other. If the multiplier is odd, add the multiplicand to the total.
It's really doing the same thing as "long multiplication" like you're taught in primary school but in binary so when you add a 0 to the right for the higher order digits you're doubling, not multiplying by ten. If you write code to do it you'd shift the multiplier first then consider whether or not to add by testing the Carry flag, or "Link bit" if like the author of the book I read you're demonstrating it on a PDP8 ;-)
But let's have a worked example, picking two numbers at random 205 * 707, use the smaller as the multiplier:
205, 707 odd, add 707 to total
102, 1414 even, disregard
51, 2828 odd, add 2828 to the total
25, 5656 odd, add 5656 to the total
12, 11312 even, disregard
6, 22624 even, disregard
3, 45248 odd, add 45248 to the total
1, 90496 odd, add 90496 to the total
--------------------------------------
144935
Yeah, that shift-and-add algorithm is sometimes used on microcontrollers, either in software if there's no hardware multiplier, or in hardware if you want the bare minimum in acceleration at a tiny cost in area.
Adds are not really considered negligible; the article is just sloppy. (Some shifts might be negligible in some models because a fixed shift requires no logic gates.) The cost of the adds in Karatsuba is significant both theoretically and in practice, and determines the cutoff where Karatsuba is useful. But the exponent in O(n^(log_2 3)) is dominated by the recursive multiplications; the adds only affect the leading constant hidden in the O().
Russian Peasant Multiplication (https://www.embeddedrelated.com/showarticle/760.php)
Mirror in Yahoo News:
https://tech.yahoo.com/science/articles/mathematicians-still...
We're going to need a bigger abacus.
Archive link. Fuck paywalls. Surprised nobody posted yet. https://archive.ph/k3dnD
Not everyone is seeing a paywall. Not sure if it’s the adblocker or other factors.
We can likely use different number representations for faster results. E.g. numbers in the form of coefficients to prime factors can be multipled at O(n) time, right?
True but addition becomes a lot less efficient in this representation :)
To make both addition and multiplication O(n), you can store numbers as their residues modulo a bunch of different primes and appeal to the Chinese Remainder Theorem. However, then size comparison becomes difficult.
Residue number systems are really neat! They're sometimes used in crypto implementations, but there you're doing modular multiplication and in most cases the modular reduction then becomes costly, so it's not a free lunch. (Except in RSA and a few other cases. RSA-CRT gets you a "free" ~4x performance boost except it's more brittle to mistakes and side-channel / fault attacks.)
There's also NTT / Fourier multiplication as an option, for big integers or polynomials or modular arithmetic.
I think the problem comes when you do a multiplication and you need more primes for uniqueness.
I think you would probably just pick enough primes at the start to handle numbers up to the number of bits you need. If we stick with primes that fit in 32-bit unsigned integers, then using the largest k such primes covers numbers up this many bits or decimal digits:
k bits digits
10 319 96
20 639 192
30 959 288
40 1279 385
50 1599 481
75 2399 722
100 3199 963
150 4799 1444
250 7999 2408
500 15999 4816
1000 31999 9632
k bits digits
10 159 48
20 319 96
30 479 144
40 639 192
50 799 240
75 1199 361
100 1598 481
150 2397 721
250 3991 1201
500 7967 2398
1000 15868 4776
k bits digits
10 78 23
20 152 45
30 220 66
40 281 84
50 327 98
54 334 100This might work really well in practice idk, but I think it's not allowed by big O to pick a maximum supported size. Otherwise you could just make a lookup table. Your algorithm must be ready for anything.
The complexity is obviously nlogn - it's just hard to prove (this comment is only somewhat serious)
If you do it in binary you only need addition.
12 × 34 = 0xC x 0x22 = 1100 x 100010
Only two 1's!
1100 add 5 zeroes + 1100 add one zero = 110011000 = 408
ta-daa!
back in 2008 I had actually a ton of success using Strassen turboplicattion kernels with some extra custom CUDA lora (its more for Cholay) than anything else but it worked. Obviously I was forbidden to use it due to interest of shana.
Karatsuba in my understanding only becomes advantageous for very large numbers relative to human scale. Mathematically it is interesting but in engineering terms the overhead usually is not worth it for practical applications. There is a fundamental trade off between factor size and product precision. If you can accept lower precision then floating point works well for large in human scale numbers.
The article discusses how Python uses it for numbers above ~2100 bits for that reason. That’s way beyond your regular floating-point type in terms of precision.
I already know about fast multiplication algorithms, but it seems there's still no proof that a faster algorithm absolutely cannot exist. In other words, we don't know where the limit is yet.
If that gets proven, would programming multiplication algorithms become faster? I'm curious
The O(n log n) algorithm is galactic (only becomes more efficient when multiplying massive numbers)
So for numbers we normally work with, no. Maybe with cryptographic operations though.
Even crypto isn't that large: 2^4096 is kind of the norm here.
Some mathematical researchers are working in the million, billion, or even trillion-bit range.
Matrix multiplication is constantly getting improved but these methods aren’t improvements on practical implementation.
Paywall.
[dead]
How is it measured? A lookup table takes 1 step to find the answer of a multiplication.
An algorithm is a finite sequence of instructions, and so can't include an infinite table. More generally, https://en.wikipedia.org/wiki/Effective_method
I disagree. 1 step is a finite sequence of instructions.
The lookup table is part of the algorithm, and is not finite.
In general any problem can be solved in 1 step with a lookup table, so here you go P=NP solved.
It doesn't have to be part of the algorithm. It all depends on how you measure it.
Well, yes, which is why they don't measure it your way, as it doesn't lead to discovering anything interesting about computation to have a shortcut like that. Or if they do it's part of a larger analysis, called an oracle machine.
I am against statements like:
A: "X people don't know how to do Y"
B: "Why not do Z?"
A: "Z is too easy and boring so they actually added more restrictions to how you are allowed to do Y so that solution doesn't count"
Maybe I'm not communicating the point clearly. In order to use a table to do the whole multiplication it has to be much larger than the largest number you would want to multiply with it. A lot of the analysis of algorithms, especially multiplication as discussed here, is about astronomically large numbers, so you don't want the existence of an even more astronomically large table as a prerequisite.
We are talking about mathematicians. We are not talking about computer engineers trying to implement a physical GPU that can multiply matrices as fast as possible. You don't have to physically build the astronomically look up table by hand. One can simply state that it exists. Proofs do not have to be clean to be a proof. You might not "want" something in a proof, but if it works then it works.