Why is my implementation of Atkin sieve is slower than Eratosthenes? [closed]

Closed. This question is not reproducible or was caused by typos. It is not currently accepting answers.

---

This question was caused by a typo or a problem that can no longer be reproduce开发者_如何学Cd. While similar questions may be on-topic here, this one was resolved in a way less likely to help future readers.

Closed last year.

Improve this question

I'm doing problems from Project Euler in Ruby and implemented Atkin's sieve for finding prime numbers but it runs slower than sieve of Eratosthenes. What is the problem?

def atkin_sieve(n)
  primes = [2,3,5]

  sieve = Array.new(n+1, false)

  y_upper = n-4 > 0 ? Math.sqrt(n-4).truncate : 1
  for x in (1..Math.sqrt(n/4).truncate) 
    for y in (1..y_upper)
      k = 4*x**2 + y**2
      sieve[k] = !sieve[k] if k%12 == 1 or k%12 == 5
    end
  end

  y_upper = n-3 > 0 ? Math.sqrt(n-3).truncate : 1
  for x in (1..Math.sqrt(n/3).truncate)
    for y in (1..y_upper)
      k = 3*x**2 + y**2
      sieve[k] = !sieve[k] if k%12 == 7
    end
  end

  for x in (1..Math.sqrt(n).truncate)
    for y in (1..x)
      k = 3*x**2 - y**2
      if k < n and k%12 == 11
    sieve[k] = !sieve[k]
      end
    end
  end

  for j in (5...n)
    if sieve[j]
      prime = true
      for i in (0...primes.length)
        if j % (primes[i]**2) == 0
          prime = false
          break
        end
      end
      primes << j if prime
    end
  end
  primes
end

def erato_sieve(n)
  primes = []
  for i in (2..n)
    if primes.all?{|x| i % x != 0}
      primes << i
    end
  end
  primes
end

---

As Wikipedia says, "The modern sieve of Atkin is more complicated, but faster when properly optimized" (my emphasis).

The first obvious place to save some time in the first set of loops would be to stop iterating over y when 4*x**2 + y**2 is greater than n. For example, if n is 1,000,000 and x is 450, then you should stop iterating when y is greater than 435 (instead of continuing to 999 as you do at the moment). So you could rewrite the first loop as:

for x in (1..Math.sqrt(n/4).truncate)
    X = 4 * x ** 2
    for y in (1..Math.sqrt(n - X).truncate)
        k = X + y ** 2
        sieve[k] = !sieve[k] if k%12 == 1 or k%12 == 5
    end
end

(This also avoids re-computing 4*x**2 each time round the loop, though that is probably a very small improvement, if any.)

Similar remarks apply, of course, to the other loops over y.

---

A second place where you could speed things up is in the strategy for looping over y. You loop over all values of y in the range, and then check to see which ones lead to values of k with the correct remainders modulo 12. Instead, you could just loop over the right values of y only, and avoid testing the remainders altogether.

If 4*x**2 is 4 modulo 12, then y**2 must be 1 or 9 modulo 12, and so y must be 1, 3, 5, 7, or 11 modulo 12. If 4*x**2 is 8 modulo 12, then y**2 must be 5 or 9 modulo 12, so y must be 3 or 9 modulo 12. And finally, if 4*x**2 is 0 modulo 12, then y**2 must be 1 or 5 modulo 12, so y must be 1, 5, 7, 9, or 11 modulo 12.

---

I also note that your sieve of Eratosthenes is doing useless work by testing divisibility by all primes below i. You can halt the iteration once you've test for divisibility by all primes less than or equal to the square root of i.

---

It would help a lot if you actually implemented the Sieve of Eratosthenes properly in the first place.

The critical feature of that sieve is that you only do one operation per time a prime divides a number. By contrast you are doing work for every prime less than the number. The difference is subtle, but the performance implications are huge.

Here is the actual sieve that you failed to implement:

def eratosthenes_primes(n)
    primes = []
    could_be_prime = (0..n).map{|i| true}
    could_be_prime[0] = false
    could_be_prime[1] = false
    i = 0
    while i*i <= n
        if could_be_prime[i]
            j = i*i
            while j <= n
                could_be_prime[j] = false
                j += i
            end
        end
        i += 1
    end
    return (2..n).find_all{|i| could_be_prime[i]}
end

Compare this with your code for finding all of the primes up to 50,000. Also note that this can easily be sped up by a factor of 2 by special casing the logic for even numbers. With that tweak, this algorithm should be fast enough for every Project Euler problem that needs you to compute a lot of primes.

---

@Gareth mentions some redundant calculations regarding 4x^2+y^2. Both here and in other places where you have calculations within a loop, you can make use of calculations you've already performed and reduce this to simple addition.

Rather than X=4 * x ** 2, you could rely on the fact that X already has the value of 4 * (x-1) ** 2. Since 4x^2 = 4(x-1)^2 + 4(2x - 1), all you need to do is add 8 * x - 4 to X. You can use this same trick for k, and the other places where you have repeated calculations (like 3x^2 + y^2).