Project Euler #58: Primes along the diagonals of a square spiral I am working on Project Euler Problem 58:


  Starting with 1 and spiralling anticlockwise in the following way, a
  square spiral with side length 7 is formed.

37 36 35 34 33 32 31
38 17 16 15 14 13 30
39 18  5  4  3 12 29
40 19  6  1  2 11 28
41 20  7  8  9 10 27
42 21 22 23 24 25 26
43 44 45 46 47 48 49

  
  It is interesting to note that the odd squares lie along the bottom
  right diagonal, but what is more interesting is that 8 out of the 13
  numbers lying along both diagonals are prime; that is, a ratio of 8/13
  ≈ 62%.
  
  If one complete new layer is wrapped around the spiral above, a square
  spiral with side length 9 will be formed. If this process is
  continued, what is the side length of the square spiral for which the
  ratio of primes along both diagonals first falls below 10%?


My code gives the correct answer, but it takes a very long time to do it, even while using Pypy, a JIT compiler.

I would like to know some good ways to increase my efficiency in this problem.

import random
def m_r(n):
    d = n - 1
    s = 0
    while d % 2 == 0:
        d &gt;&gt;= 1
        s += 1
    for repeat in range(20):
        a = 0
        while a == 0:
            a = random.randrange(n)
        if not miller_rabin_pass(a, s, d, n):
            return False
    return True

def miller_rabin_pass(a, s, d, n):
    a_to_power = pow(a, d, n)
    if a_to_power == 1:
        return True
    for i in range(s-1):
        if a_to_power == n - 1:
            return True
        a_to_power = (a_to_power * a_to_power) % n
    return a_to_power == n - 1

def eraSieve(n):
    sieve=[True]*(n+1)
    sieve[:2] = [False, False]
    sqrt = int(n**.5)+1
    for x in xrange(2, sqrt):
        if sieve[x]:
            sieve[2*x::x]=[False]*(n/x-1)
    return sieve

def diagonalNum(n): # n is the number of row
    increment = 2
    getDiaNum = 1
    limit = (n-1)/2
    sideLen = 1.0
    c = 0.0
    for i in xrange(1, limit + 1):
        count = 1
        while count &lt;= 4:
            getDiaNum += increment
            if m_r(getDiaNum):
                c += 1
            count += 1
        increment += 2
        sideLen += 4
    return c / sideLen



for x in xrange(1, 100000):
    ratio = diagonalNum(x)
    print x
    if ratio &lt; 0.10:
        break
print x
A bit of math

It is interesting to see that the corners of a square of size \$s\$ have values: \$s^2-3s+3\$, \$s^2-2s+2\$, \$s^2-s+1\$, \$s^2\$. (Also \$s^2\$ is not a prime ;-))

To be continued later.