Extended Euclidean Algorithm

I’ve started a cryp­tog­ra­phy class at DePaul and dur­ing the first lec­ture we first reviewed the Extended Euclid­ean Algo­rithm.  It pro­vides a method, using only a table, pen­cil, and paper, to eas­ily find not just the gcd of two inte­gers, but the s and t val­ues as well. Aca­d­e­m­i­cally: Let a and b be two inte­gers, at least one of which is nonzero, and let d = gcd(a,b). There exists some val­ues s,t such that d = sa + tb.

When it was cov­ered in class, the table didn’t have the k col­umn as it does in Wikipedia so I didn’t quite grasp it the first time through.  When I’m unsure of some­thing I attempt to write imple­ment it in code. So below is my quick Python imple­men­ta­tion of the Extended Euclid­ean Algo­rithm for find­ing d, s, and t, given a and b in d = gcd(a,b) = sa + tb . Not quite sure what we’ll use it for, but I have no doubt that I will learn that soon enough!

# Uses the Extended Euclid­ean Algo­rithm (Table Method) to find the gcd, s, and t
def gcd(a,b,dis­play=False):
    # Setup the first two rows of the table
    s2 = 1; t2 = 0; d2 = a; k2 = 1;
    s1 = 0; t1 = 1; d1 = b; k1 = d2/d1;
    # Com­pute Rows until the GCD is found
    while True:
        # Cal­cu­late a Table Row
        s0 = s2 — k1*s1
        t0 = t2 — k1*t1
        d0 = d2 — k1*d1
        k0 = d1/d0
        # Shift all the Row Vari­ables
        d2 = d1;    d1 = d0;
        s2 = s1;    s1 = s0;
        t2 = t1;    t1 = t0;
        k2 = k1;    k1 = k0;
        # We’ve found the GCD
        if d2%d1 == 0:
    # Print the Results
    if dis­play == True:
        print “gcd(a,b)=sa+tb: gcd(“+‘a‘+”,”+‘b‘+”)=(“+‘s1‘+”)(“+‘a‘+”)+(“+‘t1‘+”)(“+‘b‘+”) = “+‘d1‘
    # Return the results
    return d1, s1, t1

Problem #53

This prob­lem over at Project Euler turned out to be a lot sim­pler than expected.  Forty-Three other prob­lems have been solved more times than this one; yet, it was dead sim­ple. Really, there is noth­ing inter­est­ing in its solu­tion, unlike Prob­lem #97 that I posted about yes­ter­day.

It’s straight com­bi­na­torics and fac­to­ri­als. So read the prob­lem page and be amazed at what a gimme this one is :)

Result: [snip]
Time: 0.321524143219

import time
def fac­to­r­ial(n):
    if n > 1:
        return n * fac­to­r­ial(n–1)
        return 1

def com­bi­na­tions(r,n):
    return fac­to­r­ial(n)/(fac­to­r­ial(r)*fac­to­r­ial(n-r))

total = 0
start = time.time()
for n in range (1,101):
    for r in range(1, n+1):
        if com­bi­na­tions(r,n) > 1000000:
            total += 1
stop = time.time()
print “Result:”, total
print “Time:”, str(stop-start)

Problem #97

Now this prob­lem was inter­est­ing! It asked for us to find the last 10 dig­its of a very large Mersenne Prime dis­cov­ered in 2004 of the form 28433×27830457+1.

I solved this prob­lem three dif­fer­ent ways. The first was by my own logic, the other two, due to some research. My thought was that cal­cu­lat­ing 2n can be done in a loop that mul­ti­plies the pre­vi­ous value by 2 on each iter­a­tion. Rais­ing it to the 7,830,457th power though would not only take for­ever, but could be hard to store in mem­ory. So I fig­ured that we really only care about the last 10 dig­its. This meant I could mul­ti­ple the pre­vi­ous value by 2 on each iter­a­tion, but store for the next iter­a­tion a con­sist­ing of only the last ten dig­its. Python’s list syn­tax of [-10:] came in mighty handy for this.

My Method
Result: 8739992577
Time: 9.89434504509

def my_method():
    print ““
    print “My Method“
    print “==============“
    start = time.time()
    for p in range(2,7830458):
        last = str(int(last)*2)
        last = last[-10:]
    result =  str(28433*int(last)+1)
    stop = time.time()
    print “Result: “+result[-10:]
    print “Time: “+str(stop-start)

It was fairly quick, sub 10 sec­onds, but not as fast as it could be.

The next method I tried was to uti­lize the bit­wise shift oper­a­tor. I saw this for a C++ solu­tion in the forums after solv­ing it and won­dered if Python had a bit­wise shift. Turns out it does. Hav­ing a base of 2 for the pow­ers is what makes this pos­si­ble. Rather than mul­ti­ply­ing a num­ber by 2, we can shift it’s binary value to the left by one posi­tion for the same affect. CPU’s are typ­i­cally faster at bit shift­ing than they are at mul­ti­ply­ing. This solu­tion took less than 2 sec­onds! That’s quite an improvement.

Bit­wise Method
Result: 8739992577
Time: 1.6119530201

def bit­wise():
    print ““
    print “Bit­wise Method“
    print “==============“
    start = time.time()
    ans = 1
    for n in range(1,7830458):
        ans = ans « 1
        ans = ans % 10000000000
    ans = ans * 28433
    ans = ans % 10000000000
    ans = ans + 1
    stop = time.time()
    print “Result: “+str(ans)
    print “Time: “+str(stop-start)

I’m not sure how, but I ended up read­ing this arti­cle on OSIX about quickly cal­cu­lat­ing inte­ger pow­ers. It pro­posed a func­tion in C# that used recur­sion to cal­cu­late a power in O(logn) time rather than the O(n) time that the first two meth­ods imple­mented. I trans­lated the code to Python, and sure enough, it solves this prob­lem in less than half of a second!

Power Method
Result: 8739992577
Time: 0.471854925156

def get_power(a,n): # a^n
    if n==0:
        return 1
    if a==0:
        return 0
    if n%2==0:
        return get_power(a*a, n/2)
        return a*get_power(a*a, n/2)
    return 0

def power():
    print ““
    print “Power Method“
    print “==============“
    start = time.time()
    value = get_power(2,7830457)
    ans = 28433*value+1
    ans = ans % 10000000000
    stop = time.time()
    print “Result: “+str(ans)
    print “Time: “+str(stop-start)

Problem #42

Another prob­lem solved! This prob­lem intro­duced me to some util­ity func­tions in python that I did not pre­vi­ously know about. It involved the gen­er­a­tion of two lists, check­ing if a list con­tains a spe­cific value, read­ing a file, con­vert­ing a char­ac­ter into it’s ASCII num­ber, and keep­ing track of time inter­nal to the program.

Check­ing a List for Items
I searched and came across this expla­na­tion of the in oper­a­tor. Some­how I either missed it or haven’t had enough caf­feine today. I use in all the time with loop­ing con­structs like for line in file but never real­ized that it pro­duced a Boolean output.

ASCII Val­ues of Char­ac­ters
This prob­lem required the abil­ity to assign a numer­i­cal value to each let­ter. To do this I used a stan­dard func­tion named ord.

ord(“A”) # Returns 65 — the base10 ASCII Value

Con­vert­ing it to hex is sim­ple with the func­tion named hex

hex(ord(“A”)) # Returns ‘0x41’

Below is the code:


import time

# Glob­als
tri_numbers = []
total = 0

def calc_value(word):
    value = 0
    for chr in list(word):
        value += ord(chr)64
    return value

def populate_triangle_numbers(max):
    for n in range(1,max+1):

start = time.time()

file = open(“words.txt”)
words = file.read­line().split(”,”)

for word in words:
    word = word.replace(\”, “”)
    value = calc_value(word)
    if value in tri_numbers:
        total += 1
stop = time.time()

print “Count: ” + str(total)
print “Run­time: ” + str(stop — start)

Problem #40

Solved another Project Euler prob­lem today.  From the site:

An irra­tional dec­i­mal frac­tion is cre­ated by con­cate­nat­ing the pos­i­tive integers:


It can be seen that the 12^(th) digit of the frac­tional part is 1.

If d_(n) rep­re­sents the n^(th) digit of the frac­tional part, find the value of the fol­low­ing expression.

d_(1) × d_(10) × d_(100) × d_(1000) × d_(10000) × d_(100000) × d_(1000000

My orig­i­nal inten­tions were to fig­ure out an ele­gant for­mula for deriv­ing dn given some n. I thought about it for a few min­utes and fig­ured it’d be inter­est­ing to see if the straight brute force solu­tion was quick enough. So I quickly coded it up in python, and sure enough, it’s sub-second and there­fore eas­ily quick enough.

[email protected]:~/proj_euler$ time python prob40.py

real 0m0.906s
user 0m0.724s
sys 0m0.044s


list =””

for n in range(1,1000001):
    list += str(n)

answers = [list[0], list[9], list[99], list[999], list[9999], list[99999], list[999999]]
prod = 1
for n in answers:
    prod = prod * int(n)

print prod

Problem #19

Another Project Euler prob­lem has been solved! This one almost made me feel guilty as I put the high level libraries of Python to use.

How many Sun­days fell on the first of the month dur­ing the twen­ti­eth cen­tury (1 Jan 1901 to 31 Dec 2000)?

My code was short and sweet:

from date­time import date­time

d = date­time(1901,1,1)

total_sundays = 0

for yr in range(1901,2001):
    d = d.replace(year=yr)
    for m in range(1,13):
        d = d.replace(month=m)
        if(d.week­day() == 6):
            total_sundays += 1

print total_sundays

Sur­pris­ingly, the code was also very fast:

[email protected]:~$ time python prob19.py
[[answer snip]]
real 0m0.086s
user 0m0.012s
sys 0m0.000s
[email protected]:~$

Problem #52

Project Euler Prob­lem #52 was inter­est­ing and easy to solve.  The prob­lem states:

Find the small­est pos­i­tive inte­ger, x, such that 2x, 3x, 4x, 5x, and 6x, con­tain the same digits.

My method of attack was to incre­ment through the pos­i­tive inte­gers cal­cu­lat­ing the above prod­ucts for each.  The prod­ucts were then con­verted into a string, bro­ken into a list of char­ac­ters, and finally sorted. The sorted ver­sion of the lists for each prod­uct were then com­pared to see if they were all the same. Python code is below.


import sys

for x in range(126000,250000):
    a2 = list(str(x*2)); a2.sort()
    a3 = list(str(x*3)); a3.sort()
    a4 = list(str(x*4)); a4.sort()
    a5 = list(str(x*5)); a5.sort()
    a6 = list(str(x*6)); a6.sort()

    if a2 == a3 == a4 == a5 == a6:
        print “True! x=”+str(x)
        print “a2 “+.join(a2)+” “+str(x*2)
        print “a3 “+.join(a3)+” “+str(x*3)
        print “a4 “+.join(a4)+” “+str(x*4)
        print “a5 “+.join(a5)+” “+str(x*5)
        print “a6 “+.join(a6)+” “+str(x*6)
        if x % 500 == 0:
print “Com­plete”