Karatsuba multiplication

Two numbers of 65,536 digits each. Long multiplication does it with 65,536² ≈ 4.3 billion digit multiplications: each digit of one with each digit of the other. Anatoly Karatsuba showed in 1960 that so many are not needed. Below is a race right in your browser: long multiplication, Karatsuba multiplication and the browser’s built-in BigInt on longer and longer numbers. The axes are logarithmic, so the growth rate shows in the slope of the lines.

Long multiplication has a slope of about two: numbers twice as long, four times as long to multiply. Karatsuba has about 1.6: twice as long, three times as long. On short numbers long multiplication is ahead: recursion, slicing and extra additions cost more than the multiplications saved. But the lines cross, and after that the gap only grows. Change the cutoff, the length from which pieces are multiplied the school way, and watch the crossing point move.

The idea

Cut each number in half: x = x₁B + x₀, y = y₁B + y₀. Then xy = x₁y₁B² + (x₁y₀ + x₀y₁)B + x₀y₀. The middle term need not be computed with two multiplications: it equals (x₀ + x₁)(y₀ + y₁) − x₁y₁ − x₀y₀, and both subtracted terms are needed anyway. Three multiplications of halves instead of four, and repeating the trick at every level gives n^{log₂ 3} ≈ n^1.585 instead of n².

On a small example: 3141 · 2718. Halves of two digits, B = 100: x₁ = 31, x₀ = 41, y₁ = 27, y₀ = 18. Three multiplications: 31 · 27 = 837, 41 · 18 = 738 and (31 + 41)(27 + 18) = 72 · 45 = 3240. The middle term is 3240 − 837 − 738 = 1665, and the answer is 837 · 10⁴ + 1665 · 10² + 738 = 8,537,238. Multiplying by B is just a shift, and is not counted.

Why this is more than saving a quarter. Each of the three multiplications of halves is done the same way, so are their halves, and the saving repeats at every level. The time obeys the equation T(n) = 3T(n/2) + cn: additions, subtractions and shifts cost linear time. By the master theorem, 3 > 2¹, the work is concentrated in the leaves of the recursion tree, and there are n^{log₂ 3} leaves. With four multiplications of halves there would be n^{log₂ 4} = n² leaves, as many as for long multiplication: cutting in half by itself gives nothing, only the saved multiplication wins. How the work is distributed across levels on real numbers is shown by the scene in the article “Divide and conquer”.

Code

Numbers are stored as lists of “digits” in base 10,000, lowest first. Inside the recursion no carries are needed: coefficients can safely be larger than the base and even temporarily negative; this is polynomial multiplication, and the carries are done once at the end. Below the cutoff (cutoff) the pieces are multiplied the school way. The program multiplies the two halves of n!, the products 1·…·n/2 and (n/2 + 1)·…·n, checks the result against the built-in long arithmetic, and counts the digit multiplications.

     n   digits   school products   Karatsuba (cutoff 1)   (cutoff 16)
   200      375              2200                    671          1737
   500     1135             19964                   4259          8436
  1000     2569            101956                  14753         31795
  2000     5736            508464                  49383        115776

On the two halves of 2000!, 5736 decimal digits together, 1434 “digits” in base 10,000, long multiplication makes 508,464 multiplications, and Karatsuba all the way down 49,383, ten times fewer. With a cutoff of 16 there are more multiplications, 115,776, but, as we will see below, this variant runs faster.

What it costs

Time is Θ(n^{log₂ 3}) ≈ Θ(n^1.585), where n is the length of the numbers: the leaves of the recursion tree hold n^{log₂ 3} digit multiplications, and the additions and shifts over all levels together are also Θ(n^{log₂ 3}), because the work per level grows one and a half times with each level down. Memory is O(n): each level keeps numbers half as long as the previous one, and in total along a branch of the recursion it comes to at most 2n digits. The benchmark: our long multiplication, our Karatsuba with cutoff 16, and the language’s built-in multiplication:

In Python our Karatsuba overtakes long multiplication already at a few hundred digits and is about four times faster at 2048 digits. Python’s built-in multiplication is orders of magnitude faster than both, being written in C, but its slope is the same, about 1.6: each doubling of the length roughly triples the time. This is Karatsuba inside CPython. In Go the built-in math/big grows the same way, while our Karatsuba in Go overtakes long multiplication later, around four thousand digits: the inner loop of long multiplication compiles to a few machine instructions, and memory allocation in the recursion costs relatively more against that background.

The cutoff: where to switch to long multiplication

Karatsuba all the way down to a single digit makes the fewest multiplications, but it is not the fastest. At the very bottom of the recursion tree there are the most nodes, and each has a function call, slicing, a memory allocation, three additions. Multiplying a piece of a few dozen digits the school way is cheaper than breaking it up three more levels. So all real implementations switch to long multiplication below some cutoff, and the cutoff is chosen by measurement:

The cutoff depends on the machine, and the Go developers described this in detail. For math/big they measured the multiplication of numbers from 200 to 400 words long in steps of 8 at each cutoff on five processors. On most the best cutoff turned out to be about 40, on one Intel Xeon server closer to 60; they chose 40, and on that Xeon it loses to the best by about 2% on average. The same data show why there is no single ideal cutoff: for a number of 320 words all cutoffs from 21 to 40 give the same result, since the recursion goes down to pieces of 40 words, and for a number of 240 words the cutoffs from 31 to 60 are the same. So the cutoff is chosen over many lengths at once. In GMP the cutoffs are tuned for each processor type separately and chosen when the library is built.

In the first scene you can move the cutoff yourself. In JavaScript the overhead of recursion is large, and cutoff 1 loses noticeably to cutoff 16 although it makes fewer multiplications: the number of operations and the time are not the same thing.

Sums or differences, equal and unequal factors

Sums or differences. CPython computes the middle term through sums, (x₁ + x₀)(y₁ + y₀). Go and GMP through differences: x₁y₀ + x₀y₁ = x₁y₁ + x₀y₀ − (x₁ − x₀)(y₁ − y₀). The sum of the halves can turn out one digit longer than a half, a difference never; but a difference has a sign, which has to be remembered separately. Differences also have an unexpected consequence: if the halves of the numbers coincide, the third multiplication becomes a multiplication by zero and costs almost nothing. How this once spoilt a measurement for the catalogue is told in the article on complexity from measurements.

Squaring. If x = y, all three multiplications of halves are squarings too, and a square is cheaper than a product: half of the pairwise digit products need not be computed. So squaring has its own cutoff: 140 digits instead of 70 in CPython, 80 words instead of 40 in Go.

Unequal factors. Suppose one number is a hundred times longer than the other. If the long one is cut in half, the high half of the short one is zero, and the trick degenerates. CPython in this case cuts the long number into pieces as long as the short one and multiplies the short one by each piece, giving calls with equal lengths, where Karatsuba works. Go does the same. And if the second factor is very short, one or two digits, Karatsuba has nothing to do at all: such a multiplication is linear anyway.

Hence a practical conclusion for the product of many numbers. Let us compute n! in two ways with one and the same browser BigInt: multiplying in turn, each time a long number by a short one, and by a product tree, where the factors are always roughly equal in length.

In turn, the time grows roughly as n²: each of the n multiplications goes over the whole long number. As a tree, the time grows noticeably more slowly: at each level of the tree equal halves are multiplied, and fast multiplication works at full strength. Python itself computes this way: math.factorial is built as divide and conquer over halves of the range.

What comes after Karatsuba

The idea generalises. Andrei Toom in 1963 and Stephen Cook in 1966 cut numbers not into two but into k parts and managed with 2k − 1 multiplications of parts instead of k²: with three parts, five multiplications instead of nine, and growth of n^{log₃ 5} ≈ n^1.465. Arnold Schönhage and Volker Strassen in 1971 multiplied numbers through the fast Fourier transform in O(n log n log log n), and in 2019 David Harvey and Joris van der Hoeven obtained O(n log n); it is generally believed that faster is impossible, but there is no proof.

In practice the methods line up by the length of the numbers. In GMP long multiplication is followed by Karatsuba, then Toom–Cook into three, four, six and a half and eight and a half parts, and for the longest numbers FFT. CPython and Go stop at Karatsuba. The same trick works for polynomials, where there are no carries and it is even simpler; that is how the listing above writes it.

History

In the autumn of 1960 Andrey Kolmogorov conjectured at a seminar on mathematical problems of cybernetics that two n-digit numbers cannot be multiplied in fewer than on the order of n² operations. Within a week the twenty-three-year-old Anatoly Karatsuba brought a method with three multiplications of halves. In 1962 a paper on the method came out in the journal Doklady Akademii Nauk under the names of Karatsuba and Yuri Ofman; Kolmogorov submitted it, and by Karatsuba’s own account the paper became known to its author only on receiving the offprints. This story is told in more detail in the basic course chapter.

Problems

The first problem is about how many digit multiplications Karatsuba with a cutoff makes, for numbers up to a billion digits long: the recursion tree cannot be gone through, but it has only a few different subproblem sizes per level. In the second you have to find where the record of a circular genome begins, and for that you will have to see polynomial multiplication in the problem.