The Greatest Common Divisor ( GCD ) of two integers is the largest integer that divides both of these integers. However, it is not restricted to only two integers, there might be many integers. For the sake of simplicity, we try to build small programmes in Python that facilitate our task of finding GCD of only two integers in this post.

There are many methods of finding the greatest common divisor of two numbers a and b, or simply gcd( a, b). Note that gcd(0, 0) is undefined. Hence, we need to bear in mind this case while writing our program for computing GCD.

Trial-Division Approach

Let’s find the greatest common divisor of 24 and 36. The positive common divisors of these numbers are 1, 2, 3, 4, 6, and 12. Therefore, gcd(24, 36) = 12. Here we have found all possible common divisors of 24 and 36, then picked the largest one among them, 12. The first algorithm that we use for our program will be based on this idea.

The first important step is to check whether both a and b are equal to 0. If so, the program will print the warning message and redirect it to the standard error stream.

In the heart of the program, we go from 1 to a and check which number divides both a and b. As the loop condition is controlled by variable a, it is vital to make sure that a > 0. If a is 0, we return b. However, what happens if the function gcd(a,b) gets negative arguments? There will be no problem, because gcd(a,b) = gcd(−a,b) = gcd(a,−b) = gcd(−a,−b).

def gcd(a,b): # if a and b are both zero, print an error and return 0 if a == 0 and b == 0: print("WARNING: gcd called with both arguments equal to zero.", file=sys.stderr) return 0 # make sure a and b are both nonnegative if a < 0: a = -a if b < 0: b = -b # if a is zero, the answer is b if a == 0: return b d = 1 t = 1 while t <= a: if a % t == 0 and b % t == 0: d = t t += 1 return d 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 def gcd ( a , b ) : # if a and b are both zero, print an error and return 0 if a == 0 and b == 0 : print ( "WARNING: gcd called with both arguments equal to zero." , file = sys . stderr ) return 0 # make sure a and b are both nonnegative if a < 0 : a = - a if b < 0 : b = - b # if a is zero, the answer is b if a == 0 : return b d = 1 t = 1 while t <= a : if a % t == 0 and b % t == 0 : d = t t += 1 return d



As you can see computing the greatest common divisor of two integers using this algorithm is inefficient. This makes our python program very slow. We have to look for a more efficient method of finding the greatest common divisor.

Euclid’s Algorithm

This algorithm has been known since ancient times. The ancient Greek mathematician Euclid left us a description of this algorithm in his great book The Elements.

Let a = bq + r, where a, b, q, and r are integers. Then gcd(a, b) = gcd(b, r).

Proof:

Suppose that d divides both a and b. Then it follows that d also divides a − bq = r. Hence, any common divisor of a and b is also a common divisor of b and r.

Likewise, suppose that d divides both b and r. Then d also divides bq + r = a. Hence, any

common divisor of b and r is also a common divisor of a and b.

Consequently, gcd(a, b) = gcd(b, r).

Going back to our program, we make some changes to its core part, using modulo operation in Python.

def gcd(a,b): # if a and b are both zero, print an error and return 0 if a == 0 and b == 0: print("WARNING: gcd called with both arguments equal to zero.", file=sys.stderr) return 0 # make sure a and b are both nonnegative if a < 0: a = -a if b < 0: b = -b if b == 0: return a if a == 0: return b return gcd(b,a % b) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 def gcd ( a , b ) : # if a and b are both zero, print an error and return 0 if a == 0 and b == 0 : print ( "WARNING: gcd called with both arguments equal to zero." , file = sys . stderr ) return 0 # make sure a and b are both nonnegative if a < 0 : a = - a if b < 0 : b = - b if b == 0 : return a if a == 0 : return b return gcd ( b , a % b )

With the help of Euclid’s algorithm, our python program runs very fast. However, we can make the program slightly more efficient. Notice that the program is based on the recursion. To compute the greatest common divisor of two large integers may result in many calls to function gcd(a, b). The iterative solution avoids this kind of overhead and gives us a slight improvement in calculating the greatest common divisor.

def gcd(a,b): # if a and b are both zero, print an error and return 0 if a == 0 and b == 0: print("WARNING: gcd called with both arguments equal to zero.", file=sys.stderr) return 0 # make sure a and b are both nonnegative if a < 0: a = -a if b < 0: b = -b while b != 0: new_a = b new_b = a % b a = new_a b = new_b return a 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 def gcd ( a , b ) : # if a and b are both zero, print an error and return 0 if a == 0 and b == 0 : print ( "WARNING: gcd called with both arguments equal to zero." , file = sys . stderr ) return 0 # make sure a and b are both nonnegative if a < 0 : a = - a if b < 0 : b = - b while b != 0 : new_a = b new_b = a % b a = new _ a b = new_b return a



When I tested three algorithms that we have talked about so far with numbers 19,847,598 and 98,539,834 I got the following results with their corresponding time to compute the greatest common divisor:

Trial-Division: gcd = 2 Time: 5.880591 sec

Euclid-Recursion: gcd = 2 Time: 0.000010 sec

Euclid-Iteration: gcd = 2 Time: 0.000005 sec

The results are impressive, are not they?

Extended Euclid’s Algorithm

gcd(a, b) can be expressed as a linear combination with integer coefficients of a and b. These coefficients are called Bézout coefficients, named after Étienne Bézout, a French mathematician of the eighteenth.

Definition:

If a and b are positive integers, then integers s and t such that gcd(a, b) = sa + tb are called Bézout coefficients of a and b. Also, the equation gcd(a, b) = sa + tb is called Bézout’s identity.

To find Bézout’s coefficient we use Extended Euclid’s algorithm.

# return (g, x, y) a*x + b*y = gcd(x, y) def egcd(a, b): if a == 0: return (b, 0, 1) else: g, x, y = egcd(b % a, a) return (g, y - (b // a) * x, x) 1 2 3 4 5 6 7 # return (g, x, y) a*x + b*y = gcd(x, y) def egcd ( a , b ) : if a == 0 : return ( b , 0 , 1 ) else : g , x , y = egcd ( b % a , a ) return ( g , y - ( b / / a ) * x , x )

I hope my post was helpful. If you like it, please, leave a comment or share with your friends. I would be very happy!

Resources:

Discrete Mathematics and Its Applications 7th Edition – Kenneth H. Rosen

C++ for Mathematics – Edward Scheinerman

Algorithm Implementation/Mathematics/Extended Euclidean algorithm – Wikipedia

Like this: Like Loading...