Karatsuba algorithm

Karatsuba multiplication, a technique for quickly multiplying large numbres, was discovered by Karatsuba and Ofman in 1962. Its time complexity is Θ$(n^{\log _{2}3})$. This makes it faster than the classical Θ(n²) algorithm, although for small numbers it is not as fast.

Example

To multiply two n-digit numbers x and y represented in base W, where n is even and equal to 2m (if n is odd instead, or the numbers are not of the same length, this can be corrected by adding zeros at the left end of x and/or y), we can write

x = x₁ W^m + x₂

y = y₁ W^m + y₂

with m-digit numbers x₁, x₂, y₁ and y₂. The product is given by

xy = x₁y₁ W^2m + (x₁y₂ + x₂y₁) W^m + x₂y₂

so we need to quickly determine the numbers x₁y₁, x₁y₂ + x₂y₁ and x₂y₂. The heart of Karatsuba's method lies in the observation that this can be done with only three rather than four multiplications:

1. compute x₁y₁, call the result A
2. compute x₂y₂, call the result B
3. compute (x₁ + x₂)(y₁ + y₂), call the result C
4. compute C - A - B; this number is equal to x₁y₂ + x₂y₁.

To compute these three products of m-digit numbers, we can employ the same trick again, effectively using recursion. Once the numbers are computed, we need to add them together, which takes about n operations.

Complexity

If T(n) denotes the time it takes to multiply two n-digit numbers with Karatsuba's method, then we can write

T(n) = 3 T(n/2) + cn + d

for some constants c and d, and this recurrence relation can be solved, giving a time complexity of Θ(n^ln(3)/ln(2)). The number ln(3)/ln(2) is approximately 1.585, so this method is significantly faster than long multiplication. Because of the overhead of recursion, Karatsuba's multiplication is slower than long multiplication for small values of n; typical implementations therefore switch to long multiplication if n is below some threshold.

Implementations differ greatly on what the crossover point is when Karatsuba becomes faster than the classical algorithm, but generally numbers beyond the standard 2⁶⁴ "long" threshhold are faster with Karatsuba. [1] Karatsuba is surpassed by the asymptotically faster Schönhage-Strassen algorithm around 2^33,000.

An Objective Caml implementation

(* Decomposition of the numbers into arrays *)
let decomposition10 a=
  let rec log10 a=
    if a<10 then 1 else 1+(log10 (a/10))
  in
  let lga=log10 a in
  let resultat=Array.create (lga) 0 in
  let rec decomp_aux a i=
    if a=0 then () else
      begin
	let b=a/10 in
	  resultat.(i)<-(a-10*b);
	  decomp_aux b (i+1)
      end
  in
    decomp_aux a 0;
    resultat;;

(* Power function, defined recursively *)
let rec ( ** ) nb pb=
  if pb=0 then 1 else nb*(nb**(pb-1));;

(* Gives both arrays the same length *)
let equilibre vecta vectb=
  let pluslong=max (Array.length vecta) (Array.length vectb) in
  let m=Array.make_matrix 2 pluslong 0 in
    for i=0 to pluslong-1 do
      begin
	try m.(0).(i)<-vecta.(i) with _->m.(0).(i)<-0;
      end;
      begin 
	try m.(1).(i)<-vectb.(i) with _->m.(1).(i)<-0;
      end
    done;
    m.(0),m.(1);;

(* Karatsuba multiplication function *)
let karatsuba a b=
  let decompa0=decomposition10 a
  and decompb0=decomposition10 b in
  let decompa,decompb=equilibre decompa0 decompb0 in
  let rec karatsuba_aux xi xj yi yj=
    if xi=xj && yi=yj then decompa.(xi)*decompb.(yi) else
      begin
	let x1y1=karatsuba_aux xi ((xi+xj)/2+1) yi ((yi+yj)/2+1)
	and x2y2=karatsuba_aux ((xi+xj)/2) xj ((yi+yj)/2) yj
	and x1y2=karatsuba_aux xi ((xi+xj)/2+1) ((yi+yj)/2) yj
	and x2y1=karatsuba_aux ((xi+xj)/2) xj yi ((yi+yj)/2+1)
	and m2=(xi-xj)/2+1 in
	  x1y1*(10**(2*m2))+(x1y2+x2y1)*(10**m2)+x2y2
      end
  in
    karatsuba_aux (Array.length decompa -1) 0 (Array.length decompa-1) 0;;

References

• A. Karatsuba and Yu Ofman, "Multiplication of Many-Digital Numbers by Automatic Computers." Doklady Akad. Nauk SSSR Vol. 145 (1962), pp. 293–294. Translation in Physics-Doklady 7 (1963), pp. 595–596.
• Karatsuba Multiplication on MathWorld
• Bernstein, D. J., "Multidigit multiplication for mathematicians". Covers Karatsuba and many other multiplication algorithms.