4843.sagews

def distance(p,q):
    p = vector(p)
    q = vector(q)
    d = norm(p-q)
    return d
    
def checkGaussianPrimes(x,y):
    if x!= 0 and y != 0:
        d = x**2 +y**2
        if d in Primes():
            return True
    else:
        p = int(sqrt(x**2 + y**2))
        if mod((p-3),4) == 0 and p in Primes():
            return True
            
def gaussianPrimes(num):
    points = []
    for x in range(0,num+1):
        for y in range(0,num+1):
            if x>=y:
                if checkGaussianPrimes(x,y) == True:
                    points.append((x,y))      

    return points       


def stepList(stepSize):
    Slist = []    
    for x in range(-stepSize,stepSize +1 ):
        for y in range(-stepSize,stepSize +1 ):
            if x ==0 and y==0:
                pass   
            else:
                if (mod(x,2)==0 and mod(y,2)==0) or (mod(x,2)!=0 and mod(y,2)!=0):
                    d = distance((0,0),(x,y))
                    if d <=(stepSize):
                        Slist.append((x,y))
    
    return Slist
    
def primeList(point, stepList):
    Plist = []
    for apoint in stepList:
        x = point[0]+apoint[0]
        y = point[1]+apoint[1]
        if checkGaussianPrimes(x,y) == True:
            Plist.append((x,y)) 
    return Plist      
    
            
def buildPath(start, stepList):
    Glist = []
    laststage = start
    counter = 1
    currentstage = laststage
    nextstage = []
    while currentstage: 
        for point in currentstage:
            newprime = primeList(point,stepList)
            for p in newprime:
                if p not in currentstage and p not in laststage and p not in nextstage:
                    counter +=1
                    nextstage.append(p)
                    Glist.append(p)
        laststage = currentstage[:]
        currentstage = nextstage[:]
        del nextstage[:]
        nextstage = nextstage[:]
    return Glist              
                             
def findMax(alist):
    max = 0.0
    maxPoint = None
    for point in alist:
        d = distance((0,0),point)
        if d>max:
            max = d
            maxPoint = point
    return maxPoint   
    
def gaussianPlot(alist):
    g = Graphics()
    for apoint in alist:
        p = point(apoint)
        g = g+p
    return g   
    
def gaussianPlotAll(alist):
    g = Graphics()
    for apoint in alist:
        x1 , y1 = apoint[0], apoint[1]
        p1 = point((x1,y1))
        x2 , y2 = y1 , x1
        p2 = point((x2 , y2))
        x3 , y3 = -x1 , y1
        p3 = point((x3 , y3))
        x4 , y4 = y3 , x3
        p4 = point((x4 , y4))
        x5 , y5 = x1 , -y1
        p5 = point((x5 , y5))
        x6 , y6 = y5 , x5
        p6 = point((x6 , y6))
        x7 , y7 = -x1 , -y1
        p7 = point((x7 , y7))
        x8 , y8 = y7 , x7
        p8 = point((x8 , y8))
        
        g = g + p1+p2+p3+p4+p5+p6+p7+p8
    return g
     
    
def drawPath(alist,stepSize):
    g = Graphics()
    current = alist[0]
    for apoint in alist:
        for i in xrange(len(alist)):
            if distance(apoint,alist[i]) <= stepSize:
                x1 , y1 = apoint[0], apoint[1]
                p1 , q1 = alist[i][0] , alist[i][1]
                l1 = line([(x1,y1),(p1,q1)])
                '''l2 = line([(-x1,y1), (-p1,q1)])
                l3 = line([(-x1,-y1), (-p1,-q1)])
                l4 = line([(x1,-y1), (p1,-q1)])
                l5 = line([(y1,x1),(q1,p1)])
                l6 = line([(-y1,x1),(-q1,p1)])
                l7 = line([(-y1,-x1),(-q1,-p1)])
                l8 = line([(y1,-x1),(q1,-p1)])'''

                g = g + l1 

    return g
    
def main():
    g = Graphics()
    t = cputime()
    stepSize = sqrt(4)
    steps = stepList(stepSize)
    path = buildPath([(-42,17)],steps)
    graph = drawPath(path, stepSize)
    print " done graphing"
    plot = gaussianPlotAll(gaussianPrimes(5))
    g =  plot
    g.show()
    print "Time Elapsed " + str(cputime(t))



main()