The greatest common denominator algorithm is often called Euclid's Algorithm. It appears in Euclid's Elements [L1 ] (Book 7, Propositions 1 and 2 [L2 ]). KPV [Ugh; Keith, CL doesn't like that claim. It strikes him as analogous to calling Apache the Web. In my idiolect, Euclid's is an algorithm--certainly not the only one--for computing the GCD.]
KPV How about if I rephrase it thus: The most common and usually most efficient algorithm to compute the GCD is called Euclid's algorithm and is based on finding the remainder of the smaller number repeatedly subtracted from the larger. It appears in Euclid's Elements [L3 ] (Book 7, Propositions 1 and 2 [L4 ]). It may have been derived from an earlier algorithm due to Eudoxus around 375 BCE. Euclid's algorithm may hold the distinction of being the oldest nontrivial algorithm (the Russian peasant algorithm used for multiplication may be older). Binary GCD below gives a non Euclid algorithm for computing GCD.
(moved out of the Bag of algorithms) - RS started this with:
proc gcdA {p q} {
while {1} {
if {$p==$q} {
# Termination case
return $p
} elseif {$p>$q} {
# Swap p and q
set t $p
set p $q
set q $t
}
# Work out q%p ('cos we're looping...)
incr q [expr {-$p}]
}
} ;#RS, with recursion removed by DKFIt is simpler and faster to use % (reminder) instead of decrementing. (Consider what happens when the above tries to compute the GCD of 1000 and 10.) Furthermore the above control structures can be greatly simplified:
proc gcdB {p q} {
set p [expr {abs($p)}]
set q [expr {abs($q)}]
while {$q != 0} {
set r [expr {$p % $q}]
set p $q
set q $r
}
set p
}There is however an even shorter way of doing it:
proc gcdC {p q} {
set p [expr {abs($p)}]
set q [expr {abs($q)}]
while {$q != 0} {set q [expr {$p % [set p $q]}]}
set p
}/Lars Hellström
CL wants to flatten that even more:
proc gcdE {p q} {
while {$q != 0} {set q [expr {abs($p) % [set p [expr abs($q)]]}]}
set p
} KPV how about
proc gcdH {p q} {
while {$q != 0} {set q [expr {$p % [set p $q]}]}
expr {abs($p)}
}Cameron Laird wrote in the comp.lang.tcl newsgroup: Simpler:
proc gcdD {a b} {
if {$b==0} {
return $a
}
gcdD $b [expr {$a % $b}]
} proc gcdE {a b} {expr {$b==0? $a: [gcdE $b [expr {$a%$b}]]}}You can also do it iteratively, with arbitrary-precision integers, ...
Tom Poindexter added: ... e.g., by using Mpexpr, it's a one-liner:
% package require Mpexpr 1.0 % mpformat %r [mpexpr 60/100.0] 3/5
Mpexpr: ftp://ftp.procplace.com/pub/tcl/sorted/math/Mpexpr-1.0/
Extended Euclid's Algorithm kpv
It is easy to rewrite the algorithm to include additional useful information. Specifically, the integer coefficients x and y such that
d = gcd(p,q) = px + qy
proc Xgcd {p q} {
if {$q == 0} {
return [list $p 1 0]
}
foreach {d x y} [Xgcd $q [expr {$p % $q}]] break
return [list $d $y [expr {$x - $y * int($p/$q)}]]
}
On a theoretical note, the analysis of Euclid's algorithm brings up a surprising connection Fibonacci numbers. It turns out the worst case input for the algorithm are two consecutive Fibonacci numbers.
Binary GCD - a non Euclid algorithm GCD KPV
GCD(p,q) = if p and q are both even => gcd(p,q) = 2gcd(p/2,q/2)
if p is odd and q is even => gcd(p,q) = gcd(p,q/2)
if p is even and q is odd => gcd(p,q) = gcd(p/2,q)
if p and q are both odd => gcd(p,q) = gcd((p-q)/2,q)This method avoids the modulus computation found in the other methods, and only uses parity tests, subtractions and shifts. On some computers, this may be faster.
Preformance - I (BBH) changed the names of the above to be unique and ran this little performance comparison:
set nums {10 11 15 17 35 95 160 1000 1111 1113}
puts " each block contains the time for \[gcdX \$p \$q]/\[gcdX \$q \$p] for 5000 iterations"
puts " |-----------------------------------------------------------|"
puts " | P & Q | gcdA | gcdB | gcdC | gcdD |"
set i 0
foreach p $nums {
incr i
foreach q [lrange $nums $i end] {
puts " |-----------+-----------+-----------+-----------+-----------|"
puts -nonewline [format " | %4d %4d |" $p $q]
foreach f {A B C D} {
# time for p q AND for Q P to see if it makes any difference
puts -nonewline [format " %4d/%4d |" [lindex [time {gcd$f $p $q} 5000] 0] [lindex [time {gcd$f $q $p} 5000] 0]]
}
puts ""
update
}
}
puts " |-----------------------------------------------------------|"and the results on my PIII running Win2K tcl 8.3.4 on 23 Jan 2002 were
|-----------------------------------------------------------| | P & Q | gcdA | gcdB | gcdC | gcdD | |-----------+-----------+-----------+-----------+-----------| | 10 11 | 52/ 56 | 28/ 22 | 24/ 20 | 28/ 22 | |-----------+-----------+-----------+-----------+-----------| | 10 15 | 26/ 28 | 26/ 24 | 26/ 20 | 28/ 24 | |-----------+-----------+-----------+-----------+-----------| | 10 17 | 44/ 46 | 32/ 28 | 28/ 26 | 36/ 32 | |-----------+-----------+-----------+-----------+-----------| | 10 35 | 34/ 36 | 26/ 22 | 24/ 20 | 28/ 24 | |-----------+-----------+-----------+-----------+-----------| | 10 95 | 54/ 54 | 26/ 22 | 24/ 20 | 28/ 22 | |-----------+-----------+-----------+-----------+-----------| | 10 160 | 68/ 68 | 22/ 20 | 20/ 20 | 22/ 20 | |-----------+-----------+-----------+-----------+-----------| | 10 1000 | 341/ 350 | 24/ 20 | 22/ 18 | 22/ 20 | |-----------+-----------+-----------+-----------+-----------| | 10 1111 | 415/ 437 | 24/ 24 | 22/ 20 | 28/ 22 | |-----------+-----------+-----------+-----------+-----------| | 10 1113 | 401/ 403 | 28/ 28 | 24/ 22 | 32/ 28 | |-----------+-----------+-----------+-----------+-----------| | 11 15 | 44/ 46 | 30/ 30 | 28/ 24 | 36/ 30 | |-----------+-----------+-----------+-----------+-----------| | 11 17 | 46/ 50 | 32/ 28 | 28/ 24 | 36/ 34 | |-----------+-----------+-----------+-----------+-----------| | 11 35 | 54/ 52 | 30/ 26 | 24/ 24 | 30/ 28 | |-----------+-----------+-----------+-----------+-----------| | 11 95 | 70/ 70 | 34/ 32 | 28/ 28 | 40/ 36 | |-----------+-----------+-----------+-----------+-----------| | 11 160 | 90/ 92 | 32/ 28 | 26/ 24 | 34/ 32 | |-----------+-----------+-----------+-----------+-----------| | 11 1000 | 348/ 351 | 28/ 26 | 24/ 22 | 32/ 32 | |-----------+-----------+-----------+-----------+-----------| | 11 1111 | 364/ 365 | 22/ 20 | 20/ 18 | 24/ 18 | |-----------+-----------+-----------+-----------+-----------| | 11 1113 | 372/ 373 | 30/ 24 | 28/ 22 | 36/ 30 | |-----------+-----------+-----------+-----------+-----------| | 15 17 | 54/ 58 | 28/ 28 | 24/ 22 | 34/ 28 | |-----------+-----------+-----------+-----------+-----------| | 15 35 | 34/ 34 | 26/ 24 | 22/ 24 | 34/ 32 | |-----------+-----------+-----------+-----------+-----------| | 15 95 | 54/ 48 | 26/ 24 | 22/ 20 | 28/ 24 | |-----------+-----------+-----------+-----------+-----------| | 15 160 | 62/ 62 | 30/ 26 | 24/ 22 | 32/ 28 | |-----------+-----------+-----------+-----------+-----------| | 15 1000 | 246/ 244 | 30/ 24 | 26/ 22 | 32/ 26 | |-----------+-----------+-----------+-----------+-----------| | 15 1111 | 306/ 310 | 26/ 24 | 24/ 20 | 26/ 24 | |-----------+-----------+-----------+-----------+-----------| | 15 1113 | 300/ 294 | 30/ 22 | 22/ 22 | 26/ 24 | |-----------+-----------+-----------+-----------+-----------| | 17 35 | 80/ 82 | 24/ 24 | 22/ 20 | 28/ 22 | |-----------+-----------+-----------+-----------+-----------| | 17 95 | 64/ 64 | 38/ 32 | 30/ 28 | 40/ 36 | |-----------+-----------+-----------+-----------+-----------| | 17 160 | 72/ 76 | 32/ 28 | 28/ 24 | 36/ 30 | |-----------+-----------+-----------+-----------+-----------| | 17 1000 | 238/ 238 | 34/ 32 | 28/ 28 | 40/ 36 | |-----------+-----------+-----------+-----------+-----------| | 17 1111 | 258/ 260 | 32/ 28 | 26/ 26 | 36/ 30 | |-----------+-----------+-----------+-----------+-----------| | 17 1113 | 262/ 266 | 28/ 26 | 24/ 24 | 32/ 26 | |-----------+-----------+-----------+-----------+-----------| | 35 95 | 44/ 44 | 32/ 30 | 26/ 26 | 36/ 32 | |-----------+-----------+-----------+-----------+-----------| | 35 160 | 54/ 68 | 36/ 30 | 26/ 26 | 36/ 30 | |-----------+-----------+-----------+-----------+-----------| | 35 1000 | 128/ 130 | 32/ 28 | 26/ 26 | 34/ 30 | |-----------+-----------+-----------+-----------+-----------| | 35 1111 | 162/ 166 | 34/ 32 | 28/ 26 | 40/ 36 | |-----------+-----------+-----------+-----------+-----------| | 35 1113 | 136/ 138 | 30/ 24 | 26/ 22 | 32/ 26 | |-----------+-----------+-----------+-----------+-----------| | 95 160 | 54/ 54 | 32/ 30 | 26/ 24 | 36/ 32 | |-----------+-----------+-----------+-----------+-----------| | 95 1000 | 92/ 94 | 32/ 28 | 26/ 24 | 36/ 32 | |-----------+-----------+-----------+-----------+-----------| | 95 1111 | 100/ 102 | 42/ 42 | 34/ 32 | 54/ 50 | |-----------+-----------+-----------+-----------+-----------| | 95 1113 | 128/ 128 | 38/ 36 | 34/ 30 | 44/ 40 | |-----------+-----------+-----------+-----------+-----------| | 160 1000 | 50/ 52 | 26/ 22 | 24/ 20 | 28/ 22 | |-----------+-----------+-----------+-----------+-----------| | 160 1111 | 120/ 120 | 38/ 34 | 30/ 30 | 44/ 40 | |-----------+-----------+-----------+-----------+-----------| | 160 1113 | 144/ 148 | 36/ 34 | 28/ 28 | 42/ 36 | |-----------+-----------+-----------+-----------+-----------| | 1000 1111 | 475/ 475 | 28/ 36 | 44/ 24 | 32/ 32 | |-----------+-----------+-----------+-----------+-----------| | 1000 1113 | 104/ 104 | 44/ 42 | 36/ 34 | 54/ 48 | |-----------+-----------+-----------+-----------+-----------| | 1111 1113 | 2303/2143 | 32/ 26 | 26/ 28 | 34/ 26 | |-----------------------------------------------------------|
KBK is curious how they compare on 'bad' inputs such as 2178309 and 1346269 - the smallest pair that require 32 trips through the loop, by the way.
Category Mathematics - Arts and crafts of Tcl-Tk programming