Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I'm attempting to construct a program (in C++) that will count the prime factors of a given number for a Project Euler problem using the Sieve of Atkin, but I'm having trouble understanding a few parts of the pseudocode on the Wikipedia article.

I've understood what the algorithm is doing from an academic point of view by reading the sections before and after the pseudocode, but when it comes to implementing it, I'm running up against a brick wall. For reference, here's the pseudocode:

// arbitrary search limit
limit ← 1000000         

// initialize the sieve
is_prime(i) ← false, i ∈ [5, limit] 

// put in candidate primes: 
// integers which have an odd number of
// representations by certain quadratic forms
for (x, y) in [1, √limit] × [1, √limit]:
    n ← 4x²+y²
    if (n ≤ limit) and (n mod 12 = 1 or n mod 12 = 5):
        is_prime(n) ← ¬is_prime(n)
    n ← 3x²+y²
    if (n ≤ limit) and (n mod 12 = 7):
        is_prime(n) ← ¬is_prime(n)
    n ← 3x²-y²
    if (x > y) and (n ≤ limit) and (n mod 12 = 11):
        is_prime(n) ← ¬is_prime(n)

// eliminate composites by sieving
for n in [5, √limit]:
    if is_prime(n):
        // n is prime, omit multiples of its square; this is
        // sufficient because composites which managed to get
        // on the list cannot be square-free
        is_prime(k) ← false, k ∈ {n², 2n², 3n², ..., limit} 

print 2, 3
for n in [5, limit]:
    if is_prime(n): print n

I'm unsure as to why the square root of limit variable is being used, rather than limit itself. I'm also unsure why the loop that removes composites by sieving is running between 5 and sqrt(limit) rather than 1 and limit.

share|improve this question
    
Jacob, 20 is divisible by 2... so no need to verify 10 because 20/2=10... every number after the square of 20 that 20 can divide is a result of division of 20/n<=sqrt(20)... –  user45754 Oct 23 '12 at 14:35

1 Answer 1

up vote 5 down vote accepted

For optimization, the algorithm treats 2 and 3 as special, and since we know that 4 isn't a prime, the latter loop starts at 5.

You also don't need to go over sqrt(limit) (a critical optimization in every prime finding algorithm) because , , and will always be larger than the limit itself (A number can never be evenly divisible by a number bigger than its square root without being divisible by a smaller number).

The sieving loop should of course also at least increment 2 every iteration, since no even number over 5 can be prime.

The wiki page also tells you that the first bounds aren't optimized. For the first part, $4x^2+y^2 \leq \mathrm{limit}$, so $y \leq \sqrt{\mathrm{limit}-4x^2}$ would be a tigher bound for the inner loop (while also saving you the computational cost to check that $n\leq\mathrm{limit}$. More aggressive optimizations should be applied to get a fast implementation.

I don't see anything in your question which should prove difficult in implementation.

share|improve this answer
1  
"A number can never be evenly divisible by a number bigger than its square root)" are you suggesting that 20 isn't evenly divisible by 10? –  JSchlather Aug 22 '11 at 17:47
    
Thanks for spotting that. I hope i didn't make any more embarrassing mistakes like that. –  Mikael Öhman Aug 22 '11 at 17:55

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

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