"""
Euclid (Book IX, Prop 36.):  If M is a
Mersenne prime, then the Mth triangular
number is a perfect number.

Context:
http://www.4dsolutions.net/ocn/perfectnumbers.html

(gpl) 2007 by K. Urner, 4D Solutions
"""

# Table 5.1 The Book of Numbers, JH Conway, RK Guy.
mersennep   = [2, 3, 5, 7, 13, 17, 19,
              31, 61, 89, 107, 127, 521, 607,
              1279, 2203, 2281, 3217, 4253, 4423, 9689,
              9941, 11213, 19937, 21701, 23209, 44497,
               86243]


mersennes = (pow(2,p) - 1 for p in mersennep )

def trinumb(n):
    """
    Triangular numbers
    """
    return n*(n+1)//2

def divisors(n):
    """
    Divisors of small numbers:  trial by division

    Changed slightly since video (see web page, URL above)
    """
    thedivisors = [1]
    d = 2
    while True:
        if not n % d:  # if there's no remainder...            
            thedivisors.append(d)
            otherfactor = n//d
            if d != otherfactor:
                thedivisors.append(otherfactor)
        d += 1                
        if d*d > n:
            break
    thedivisors.sort()
    return thedivisors

def test():
    """
    loop through generated Mersennes testing
    Euclid's Prop IX:36
    """
    for M in mersennes:
        print "The Mersenne is... %s" % M
        tri = trinumb(M)
        print "The Triangular Number is... %s" % tri
        print "The divisors of %s are...\n" % tri
        thedivisors = divisors(tri)
        print thedivisors
        print "...and the sum of these divisors is"
        print "(drum roll)..."
        print sum(thedivisors)
        print
        yield None