import random

class Point:

    def __init__(self, x, y):
        self.x = x
        self.y = y

    def __str__(self):
        return "(%d, %d)" % (self.x, self.y)

def pointsEqual(p1, p2):
    return p1.x == p2.x and p1.y == p2.y

def egcd(a, b):
    if a == 0:
        return (b, 0, 1)
    else:
        g, y, x = egcd(b % a, a)
        return (g, x - (b // a) * y, y)

def modinv(a, m):
    g, x, y = egcd(a, m)
    if g != 1:
        raise Exception('modular inverse for %d does not exist mod %d' % (a, m))
    else:
        return x % m

def addH(p1, p2, n):
    # Hint: the back-door is entirely in this function
    abcd = (p1.x*p1.y*p2.x*p2.y) % n
    x3 = (abcd*modinv((p1.y*p2.y - p1.x*p2.x) % n, n)) % n
    y3 = (abcd*modinv((p1.x*p2.y + p1.y*p2.x) % n, n)) % n
    return Point(x3, y3)

def doubleH(p, n):
    return addH(p, p, n)

def mulH(m, p, n):
    r = None
    while m != 0:
        if m & 1:
            if r == None:
                r = p
            else:
                r = addH(r, p, n)
        m >>= 1
        p = doubleH(p, n)
    return r

# A 255 bit prime
n = 57896044618658097711785492504343953926634992332820282019728792003956564819949
# A random point picked to make you go "Hmm...."
g = Point(9, 37562988266855220659555181903926734170782983985553212517523748423444692139187)
print "g = ", g

# Alice's secret
sa = random.randrange(0, n)
# Bob's secret
sb = random.randrange(0, n)

print "sa:", sa
print "sb:", sb
A = mulH(sa, g, n)
B = mulH(sb, g, n)
print "A:", A
print "B:", B
# Shared master secret computed by Alice
aliceMS = mulH(sa, B, n)
# Shared master secret computed by Bob
bobMS = mulH(sb, A, n)
print "Alice:", aliceMS
print "Bob:", bobMS
assert pointsEqual(aliceMS, bobMS)