3
\$\begingroup\$

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 >>= 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 <= 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 < 0.10:
        break
print x
\$\endgroup\$
3
  • 2
    \$\begingroup\$ Each time you call diagonalNum from the main loop, it starts again from 1, repeating all the work done on the previous call (and the call before that, and so on). Why not maintain a running count? \$\endgroup\$ Sep 27, 2016 at 19:45
  • \$\begingroup\$ @GarethRees What do you mean by running count? \$\endgroup\$ Sep 27, 2016 at 20:04
  • 1
    \$\begingroup\$ @GarethRees That looks like it should be an answer. \$\endgroup\$ Sep 27, 2016 at 20:07

1 Answer 1

1
\$\begingroup\$

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 ;-))

My solution

Not an actual review but if you find this useful:

def euler58(ratio=0.1):
    """Solution for problem 58."""
    # First analysis in euler28
    # Corners are s*s-3s+3, s*s-2s+2, s*s-s+1, s*s
    # s*s is likely not to be prime
    nb_prime = 0
    for s in itertools.count(3, 2):
        for i in range(1, 4):
            if is_prime(s * s - i * s + i):
                nb_prime += 1
        if nb_prime < ratio * (2 * s - 1):
            return s
\$\endgroup\$
4
  • \$\begingroup\$ Thanks for the input! I didn't realize those equations for the points of the square. Thanks for letting me know! \$\endgroup\$ Sep 27, 2016 at 20:47
  • \$\begingroup\$ What were you going to add? You never finished \$\endgroup\$
    – Some Guy
    Feb 23, 2022 at 3:59
  • \$\begingroup\$ @SomeGuy Going point! I really can't tell what I had in mind back then. I've added my solution in case it is useful to anyone \$\endgroup\$
    – SylvainD
    Feb 23, 2022 at 8:14
  • \$\begingroup\$ Nice. Also, nn isn't just not likely to be prime, it's IMPOSSIBLE for nn to be prime. 1 isn't a prime number. \$\endgroup\$
    – Some Guy
    Feb 23, 2022 at 15:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.