Karatsuba algorithm: Difference between revisions

Revision as of 01:35, 14 June 2007

Revision as of 14:36, 14 June 2007

Line 4:

The speed of this algorithm relies partly on the ability to do [[arithmetic shift|shift]]s faster than a standard multiply, with shifts typically taking [[linear time]] on a practical computer.

We assume two ''n''-digit numbers ''x'' and ''y'', represented in [[radix|base]] ''B'', where ''B'' is a convenient value to shift with, like [[binary numeral system|two]] in mostmodern binary computers [[as of 2007|today]]. Choosing a number ''m'' less than ''n'', we can write

:''x'' = ''x''₁''B''^''m'' + ''x''₂

Line 27:

Thus ''xy'' = ''X'' ''B''^2m + ''Y'' + ''Z'' ''B''^m; we have cut down the number of multiplications from four to three.

To compute these three products of ''m''-digit numbers, we can [[recursion|recursively]] employ Karatsuba multiplication again, effectively using [[recursion]]. The shifts, additions and subtractions take linear time respectively, and can be considered negligible.

Karatsuba multiplication works best when the two operands are of similar size; thus while one can theoretically choose ''m'' arbitrarily, the best time bound would be achieved by choosing ''m'' such that the subproducts are roughly equal, that is setting ''m'' to be about ''n''/2.

Line 49:

Implementations differ greatly on what the crossover point is when Karatsuba becomes faster than the classical algorithm, but generally numbers over 2³²⁰ are faster with Karatsuba. [http://gmplib.org/manual/Karatsuba-Multiplication.html][http://ozark.hendrix.edu/~burch/proj/karat/comment1.html] [[Toom-Cook multiplication]] is a faster generalization of Karatsuba's, and is further surpassed by the fastest known, the [[Schönhage-Strassen algorithm]].

<!--

[[[ What's this doing here? ]]]

==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==