Talk:Sieve of Atkin

From Wikipedia, the free encyclopedia

I'm hoping that someone, if I don't have time to do it, will fill in details like how the sieve determines how many values for x and y fulfill the deletion question. If we can move the article more toward the easy reading of Sieve of Eratosthenes, that would be great. — 131.230.133.185 5 July 2005 19:26 (UTC)


I concure!

Contents

[edit] Class One, Class Two, Class Three

I just did some rewording on the algorithm. However, I do not know exactly about the algorithm, and am confused by the following:

Then, for each remaining number, n:
In class one, count the number of combinations of x > 0 and y > 0 that are solutions to 4x2 + y2 = n.
In class two, count the number of combinations of x > 0 and y > 0 that are solutions to 3x2 + y2 = n.
In class three, count the number of combinations of x > 0, y > 0, and x>y that are solutions to 3x2 - y2 = n.
If the count was even, delete the number from the list.

What are classes? The article says "If the count was even". If what count was even (any one of them or the sum of them)? Clarification would make the article easier to read and understand. --mets501 04:12, 4 March 2006 (UTC)

I think the "classes" are the backup lists (i.e., the ones that hold values that are not yet determined to be prime) that were just mentioned previously. I'm not sure enough about this to edit the page myself, however, as I'm not familiar with the Sieve of Atkin. In regards to the count being even, it refers to the number of combinations for each listed "class". If the number of combinations (say for class 1, of 4x^2 + y^2 = n) is even for the current n, then you get rid of the number. Any word from an expert on this would be great, however. Opblaaskrokodil 03:27, 26 March 2006 (UTC)


I would like to know if the ordere pair solutions have to be integers. For example, the non prime number 91 would be in list 2. I can find only one ordered integer pair (3, 8) to solve 3x^2 + y^2 = 91. Thus I do not remove it from the list, and it would erroneously go into the list of primes. the preceding comment is by 24.66.94.140 - 16:55, 15 April 2006 (UTC): Please sign your posts!

The solutions certainly have to be integers, otherwise there would be an infinite number of solutions for every number. As for 91, you missed (5, 4). -- EJ 22:06, 4 May 2006 (UTC)

[edit] Clean up tag removed

I removed the clean up tag because I couldn't figure out why it was put on in the first place. Anyone who wants to put the tag back on had better explain not only why but also why they couldn't fix what they thought was wrong with the page themselves in the first place. Anton Mravcek 20:48, 8 May 2006 (UTC)

I'm not sure why a cleanup tag was justified either (unless because the explanation of the algorithm leaves a bit to be desired). I did add an expandsect tag to the "proof and explanations" section; there is currently no explanation of why the business with calculating the various numbers of combinations of x and y actually helps weed out non-primes. Maybe it would be obvious to a mathematician why this is so, but it sure isn't obvious to me. In fact, I have no idea (if I did, I'd add the extra explanation myself). Zoicon5 18:06, 26 May 2006 (UTC)
I went and looked up the original paper using the external link; this was clear enough for me to be able to add further explanation of the part I hadn't understood before, so I did that (and removed the "expandsect" tag). Zoicon5 19:13, 26 May 2006 (UTC)
I've cleaned up the article quite a bit, making it simpler and easier to follow (like the Sieve of Eratosthenes article) and showing where to obtain proofs of the various assertions. I've removed quasi-optimizations in favor of understandability (the reader can choose his own optimizations or copy the C code linked at the bottom). — 131.230.133.186 18:53, 28 May 2006 (UTC)


[edit] Confused about x and y

I'm confused about the X and Y variables appearing out of nowhere in the psuedo code? I don't know enough maths to follow well, i'm using this to implement a python script that will do this, basically. I presume the X and Y should be I and J (from the for....do... loops?) someone with more math/algorithm nous?
Well, the previous pseudocode declared x and y, but someone declared it "crud" and "wrong"; it looked like this:
HighestFactor = 1000000;               // Enough. A prime would be better.
Limit = HighestFactor^2;               // Resultant search span.
array is_prime[5:Limit] initial false; // Sieve array.
biginteger x,y,n,k;                    // Must be able to hold 5*Limit: 4x^2 + y^2!
// put in candidate primes: integers which have an odd number of representations by certain quadratic forms.
for x:=1:HighestFactor do
 for y:=1:HighestFactor do
   n:=4*x^2 + y^2; if n <= Limit and n mod 12 = 1 or = 5 then     is_prime[n]:=not is_prime[n];
   n:=3*x^2 + y^2; if n <= Limit and n mod 12 = 7 then            is_prime[n]:=not is_prime[n];
   n:=3*x^2 - y^2; if n <= Limit and x > y and n mod 12 = 11 then is_prime[n]:=not is_prime[n];
 next y;
next x;
// eliminate composites by sieving
// if n is prime, omit all multiples of its square; this is sufficient because 
// composites which managed to get on the list cannot be square-free
for n:=5:HighestFactor do
  if is_prime[n] then for k:=n^2:Limit:n^2 do is_prime[k]:=false; next k;
next n;
// Present the results.
print 2, 3; for n:= 5:Limit do if is_prime[n] then print n; next n;
And notice that the variables do not appear "out of nowhere" as you complained. As for the "wrongness", I transliterated it into Pascal, producing
Program PrimesByAtkin; uses DOS, crt;
 const HighestFactor = 100;             {Enough. A prime would be better.}
 const Limit = HighestFactor*HighestFactor;             {Resultant search span.}
 var is_prime: array[5..Limit] of boolean;              {Sieve array.}
 var x,y,n,k: longint;                  {Must be able to hold 5*Limit: 4x^2 + y^2!}
 var outa: text; {Bizarre declaration of an output file.}
 BEGIN
  for n:=5 to Limit do is_prime[n]:=false;
 {Put in candidate primes: integers which have an odd number of representations by certain quadratic forms.}
  for x:=1 to HighestFactor do
   for y:=1 to HighestFactor do
    begin
     n:=4*x*x + y*y; if (n <= Limit) and ((n mod 12 = 1) or (n mod 12 = 5)) then is_prime[n]:=not is_prime[n];
     n:=3*x*x + y*y; if (n <= Limit) and (n mod 12 = 7) then                     is_prime[n]:=not is_prime[n];
     n:=3*x*x - y*y; if (n <= Limit) and (x > y) and (n mod 12 = 11) then        is_prime[n]:=not is_prime[n];
    end; {next y;}
  {next x;}
 { eliminate composites by sieving
 if n is prime, omit all multiples of its square; this is sufficient because
 composites which managed to get on the list cannot be square-free}
  for n:=5 to HighestFactor do
   if is_prime[n] then
    begin
     {for k:=n^2:Limit:n^2 do is_prime[k]:=false; next k;}
     k:=n*n;
     repeat
      is_prime[k]:=false;
      k:=k + n*n;
     until k > Limit;
    end;
 {next n;}
 {Present the results.}
  assign(outa,'PrimesA.txt'); rewrite(outa);
  writeln(outa,2);writeln(outa,3); for n:= 5 to Limit do if is_prime[n] then writeln(outa,n); {next n;}
  close(outa);
 END.
Which is not pseudocode but actual Pascal source, which compiles and works. You should have no trouble converting into python, or anything else. Incidentally, There is no need to try "if x == true then" in a C-gibberish form, "if x then" should do.

NickyMcLean 20:10, 29 November 2006 (UTC)

[edit] How is this

Perhaps the following pseudocode example would be more suiting to you McLean. I apologize for my previous mistakes in not changing x and y to i and j respectively.

Please note that in your example:

if is_prime[n] then for k:=n^2:Limit:n^2 do

is very much gibberish

limit = 1000000                         // Resultant search span. ^ being the exponentiation operation
highest_factor = sqrt(limit)            // Square root of the search span
is_prime[limit]                         // Sieve array. Initialize all entries to the equivalent of a false boolean
(n, i, j) = (0, 0, 0)                   // Must be able to hold 5*Limit!

// Put in candidate primes: Integers which have an odd number of representations by certain quadratic forms.

for i from 1 to highest_factor do
    for j from 1 to highest_factor do
       n = 4*i^2 + j^2
       if n <= limit and ((n mod 12 == 1) or (n mod 12 == 5))  then
           is_prime[n] = not is_prime[n]
       end if

       n = 3*i^2 + j^2
       if n <= limit and (n mod 12 == 7) then
           is_prime[n] = not is_prime[n]
       end if

       n = 3*i^2 - j^2
       if n <= limit and (i > j) and (n mod 12 == 11) then
           is_prime[n] = not is_prime[n]
       end if
    end for
end for

// Eliminate composites by sieving
// If n is prime, omit all multiples of its square; this is sufficient because
// composites which managed to get on the list cannot be square-free

for i from 5 to highest_factor do
    if is_prime[i] == true then
       for j from i^2 to limit step i^2 do
           is_prime[j] = false
       end for
    end if
end for

// With the exception of 2 and 3, is_prime now contains boolean values
// for the primality of integers from 5 to limit at their respective indexes in the array

Here is a Python example, so you can see how easily my pseudocode translates into the language:

   import math

   limit = 100
   highest_factor = int(math.floor(math.sqrt(limit)))
   is_prime = [False] * limit
   # i, j and n will be dynamically created

   for i in xrange(1, highest_factor):
       for j in xrange(1, highest_factor):
           n = 4*i**2+j**2
           if n <= limit and ((n % 12 == 1) or (n % 12 == 5)):
               is_prime[n] = not is_prime[n]
           
           n = 3*i**2+j**2
           if n <= limit and (n % 12 == 7):
               is_prime[n] = not is_prime[n]
               
           n = 3*i**2 - j**2
           if n <= limit and (i > j) and (n % 12 == 11):
               is_prime[n] = not is_prime[n]
               
   for i in xrange(5, highest_factor):
       if is_prime[i]:
           for j in xrange(i**2, limit, i**2):
               is_prime[j] = False

   # Print the integers
   for i in xrange(5, limit):
       if is_prime[i]:
           print i


No, I was not the person who declared it 'wrong' or whatever, I was in fact, the person who removed that comment. I just think my pseudocode is cleaner and easier to understand. I do not mean it as an insult to you.

65.33.145.252 03:52, 30 November 2006 (UTC)

(I wasn't accusing you: the other person merely has some other IP address recorded. I was peeved by the crass remark, but against the possibility of a miscode, the Pascal transliteration offered some support though the many minor changes could introduce error)
Ah, I see now that Python is indeed one of the languages that uses indentation rather than do ... end, begin ... end or {...} (which in Pascal means a comment) - I'd half-remembered that. However I'd suggest as for layout, having the three action lines one above the other facilitates comparison and the separation of the conditions from the boringly simple actions (whose expression requires many symbols just to flip a bit) which are each the same. There are some (very obscure) languages that do use a 2-dimensional layout, but they're a different world. As a general point, depending on layout is risky because many times your carefully-indented text will be simplified to flush-left by over-helpful web page handlers. Thus it can help to signify the scope of blocks by explicit markers, which is also a hint to those unfamiliar with the notion that a visual hint (layout) does in this case translate to a formal specification of block scope.
I agree that for i:=start:stop:step do is obscure, but for i:=start:stop is common (the idea has reached Matlab, for example) and reasonably similar to a common form of array bound specification. Pascal's Array[start..stop] is particularly lame. So yes, I should have explained that, or used the form with keywords, as in your Pythonish version. I'd suggest though that the Pythonish "with" keyword is unfortunate, and could better be replaced by "step" for pseudocode, even if in Python it has some special other meaning. Ah, I see that in actual Python, a peculiar form is used, without hints. I think that the syntax for i:=start:stop:step do ... next i; (with the :step optional) is completely satisfactory, and if I ruled the world... NickyMcLean 19:50, 30 November 2006 (UTC)
Perhaps then we can work together on the pseudocode, to generate a clean and readable format using both our styles that is not only satisfactory to us both, but much easier for those trying to translate it into their preferred language and either of our examples would've been on their own. I'm more than happy to compromise 'with' for 'step' or 'next'. I look forward to your thoughts.

65.33.145.252 21:58, 30 November 2006 (UTC)

Well, I'd allow for i from n^2 to Limit step n^2 do instead of for i:=n^2:Limit:n^2 do in the spirit of easy readability for those not familiar with the form, but introducing the "cup" and "cap" symbols of formal logic for "or" and "and" without explanation seems a bit obscure. NickyMcLean 21:00, 4 December 2006 (UTC)

[edit] Limits

Everyone seems locked into the concept that in exhaustive testing for primality of N, factors up to √N must be considered, and they persist in using a name called Limit for their description of the range of N to be worked through so that they repeatedly talk about √Limit in controlling matters. Since most integers do not have integer square roots this is messy at best and in terms of computational effort, wastefully requires the evaluation of square roots. A better description of the bounds on the possible factors f would be for f^2 <= N rather than f <= √N. And then if you think about the choice of possible factors (as in the division method http://mersennewiki.org/index.php/Trial_factoring) you soon see that only primes should be considered, and it is best to start with the smallest (since more numbers are divisible by 2 than by 3, etc.) working from a small array of prime numbers, and then you need recall only the highest index so far needed, that is f = Prime(i) for i up to LastIndex, the last index for which Prime(LastIndex)^2 <= N, or similar schemes.

And then you might realise that the actual requirement to check N for primality by trial division is as follows: use as factors the primes P(i) with i from 1 to Last such that N < P(Last + 1)^2. For instance, suppose that 5 is the last prime factor being considered as a possible divisor for N. On the face of it, as soon as N reaches 25, you would extend the test list from {2,3,5} to {2,3,5,7}. But this is not necessary. There is no need to test with 7 until 49 = 7*7 is reached, because until then, N = 2*7 need not be tested against 7 because 2 would already be tested for, nor need it be tried for N = 3*7, nor for N = 5*7. Only when N = 7*7 is reached would the test divisions by {2,3,5} not suffice.

In the current context, if MaxFactor is the largest factor that will be considered, then the upper bound for N is MaxFactor^2 and there is no confusion of square roots, especially of numbers that are not perfect squares. Further, if MaxFactor is a prime, say P(i), then the upper bound will be P(i + 1)^2 - 1, though P(i)^2 will do, as would (P(i) + 2)^2 since P(i) + 2 <= P(i + 1) for i > 1. Whichever, the idea is to avoid littering the description with √N and similar. NickyMcLean 21:00, 4 December 2006 (UTC)