Pregunta

Estoy haciendo un proyecto en este momento y necesito un método eficiente para calcular números primos. He usado el tamiz de Eratóstenes pero, he estado buscando y he encontrado que el tamiz de Atkin es un método más eficiente. Me ha resultado difícil encontrar una explicación (¡que he podido entender!) De este método. ¿Como funciona? El código de ejemplo (preferiblemente en C o python) sería brillante.

Editar: gracias por su ayuda, lo único que todavía no entiendo es a qué se refieren las variables x e y en el pseudocódigo. ¿Podría alguien arrojarme algo de luz sobre esto?

¿Fue útil?

Solución

La página wiki siempre es un buen lugar para comenzar, ya que explica el algoritmo en completo y proporciona pseudocódigo comentado. (N.B. Hay muchos detalles, y dado que el sitio web wiki es confiable, no lo citaré aquí.)

Para referencias en los idiomas específicos que mencionó:

Espero que eso ayude.

Otros consejos

El artículo de Wikipedia tiene una explicación:

  • " El algoritmo ignora completamente cualquier número divisible por dos, tres o cinco. Todos los números con un resto par módulo sesenta son divisibles por dos y no primos. Todos los números con módulo sesenta resto divisible por tres también son divisibles por tres y no primos. Todos los números con módulo restante sesenta divisible por cinco son divisibles por cinco y no primos. Todos estos restos son ignorados. & Quot;

Comencemos con los famosos

primes = sieve [2..] = sieve (2:[3..])   
  -- primes are sieve of list of 2,3,4... , i.e. 2 prepended to 3,4,5...
sieve (x:xs) = x : sieve [y | y <- xs, rem y x /= 0]   -- set notation
  -- sieve of list of (x prepended to xs) is x prepended to the sieve of 
  --                  list of `y`s where y is drawn from xs and y % x /= 0

Veamos cómo procede con algunos primeros pasos:

primes = sieve [2..] = sieve (2:[3..]) 
                     = 2 : sieve p2     -- list starting w/ 2, the rest is (sieve p2)
p2 = [y | y <- [3..], rem y 2 /= 0]     -- for y from 3 step 1: if y%2 /= 0: yield y

p2 no debe contener múltiplos de 2 . IOW solo contendrá números coprime con 2 & # 8212; ningún número en p2 tiene 2 como factor. Para encontrar p2 en realidad no necesitamos probar dividir entre 2 cada número en [3 ..] (eso es 3 y arriba, 3,4,5,6,7, ... ), porque podemos enumerar todos los múltiplos de 2 por adelantado:

rem y 2 /= 0  ===  not (ordElem y [2,4..])     -- "y is not one of 2,4,6,8,10,..."

ordElem es como elem (es decir, prueba de elemento ), simplemente asume que su argumento de lista es ordenado, lista creciente de números, para que pueda detectar la no presencia de forma segura, así como la presencia:

ordElem y xs = take 1 (dropWhile (< y) xs) == [y]   -- = elem y (takeWhile (<= y) xs) 

El elem ordinario no asume nada, por lo que debe inspeccionar cada elemento de su argumento de lista, por lo que no puede manejar listas infinitas. ordElem can. Entonces,

p2 = [y | y <- [3..], not (ordElem y [2,4..])]  -- abstract this as a function, diff a b =
   = diff      [3..]                 [2,4..]    --       = [y | y <- a, not (ordElem y b)]
   -- 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
   -- . 4 . 6 . 8 . 10  . 12  . 14  . 16  . 18  . 20  . 22  .
   = diff [3..] (map (2*)            [2..] )  --  y > 2, so [4,6..] is enough
   = diff [2*k+x | k <- [0..],  x <- [3,4]]   -- "for k from 0 step 1: for x in [3,4]:
          [2*k+x | k <- [0..],  x <- [  4]]   --                            yield (2*k+x)"
   = [     2*k+x | k <- [0..],  x <- [3  ]]   -- 2 = 1*2 = 2*1
   = [3,5..]                                  -- that's 3,5,7,9,11,...

p2 es solo una lista de todas las probabilidades por encima de 2 . Bueno, sabíamos eso . ¿Qué pasa con el siguiente paso?

sieve p2 = sieve [3,5..] = sieve (3:[5,7..]) 
                         = 3 : sieve p3
p3 = [y | y <- [5,7..], rem y 3 /= 0]
   = [y | y <- [5,7..], not (ordElem y [3,6..])]           -- 3,6,9,12,...
   = diff [5,7..] [6,9..]         -- but, we've already removed the multiples of 2, (!)
   -- 5 . 7 . 9 . 11 . 13 . 15 . 17 . 19 . 21 . 23 . 25 . 27 .
   -- . 6 . . 9 . . 12  . . 15 . . 18  . . 21 . . 24  . . 27 .
   = diff [5,7..] (map (3*) [3,5..])                       -- so, [9,15..] is enough
   = diff [2*k+x | k <- [0..], x <- [5]] (map (3*)
          [2*k+x | k <- [0..], x <- [    3]] )
   = diff [6*k+x | k <- [0..], x <- [5,7,9]]               -- 6 = 2*3 = 3*2
          [6*k+x | k <- [0..], x <- [    9]] 
   = [     6*k+x | k <- [0..], x <- [5,7  ]]               -- 5,7,11,13,17,19,...

Y el siguiente,

sieve p3 = sieve (5 : [6*k+x | k <- [0..], x <- [7,11]])
         = 5 : sieve p5
p5 = [y | y <-        [6*k+x | k <- [0..], x <- [7,11]], rem y 5 /= 0]
   = diff [ 6*k+x | k <- [0..], x <- [7,11]]   (map   (5*)
          [ 6*k+x | k <- [0..], x <- [                  5,       7]]) -- no mults of 2 or 3!
   = diff [30*k+x | k <- [0..], x <- [7,11,13,17,19,23,25,29,31,35]]  -- 30 = 6*5 = 5*6
          [30*k+x | k <- [0..], x <- [                 25,      35]]
   = [     30*k+x | k <- [0..], x <- [7,11,13,17,19,23,   29,31   ]]

Esta es la secuencia con la que está trabajando el tamiz de Atkin. No hay múltiplos de 2, 3 o 5 presentes en él. Atkin y Bernstein trabajan módulo 60 , es decir, duplican el rango:

p60 = [ 60*k+x | k <- [0..], x <- [1, 7,11,13,17,19,23,29,31, 37,41,43,47,49,53,59]]

A continuación, usan algunos teoremas sofisticados para conocer algunas propiedades de cada uno de esos números y tratar cada uno en consecuencia. Todo se repite (a-la la "rueda") con un período de 60 .

  • " Todos los números n con módulo sesenta restantes 1, 13, 17, 29, 37, 41, 49 o 53 (...) son primos si y solo si el número de soluciones para 4x ^ 2 + y ^ 2 = n es impar y el número es cuadrado. & Quot;

¿Qué significa eso? Si sabemos que el número de soluciones a esa ecuación es impar para tal número, entonces es primo si es libre de cuadrados. Eso significa que no tiene factores repetidos (como 49, 121, , etc.).

Atkin y Bernstein usan esto para reducir el número total de operaciones: para cada primo (desde 7 y arriba), enumeramos (y marcamos para eliminar) los múltiplos de su cuadrado , por lo que están mucho más separados que en el tamiz de Eratóstenes, por lo que hay menos de ellos en un rango determinado.

El resto de las reglas son:

  • " Todos los números n con módulo sesenta restantes 7, 19, 31 o 43 (...) son primos si y solo si el número de soluciones para 3x ^ 2 + y ^ 2 = n es impar y el número es cuadrado. "

  • " Todos los números n con módulo sesenta restantes 11, 23, 47 o 59 (...) son primos si y solo si el número de soluciones para 3x ^ 2 & # 8722; y ^ 2 = n es impar y el número es cuadrado. "

Esto se encarga de los 16 números principales en la definición de p60 .

vea también: ¿Cuál es el algoritmo más rápido para encontrar números primos?


La complejidad temporal del tamiz de Eratóstenes para producir primos de hasta N es O (N log log N) , mientras que la del tamiz de Atkin es O (N) (dejando de lado las optimizaciones log log N adicionales que no se escalan bien). Sin embargo, la sabiduría aceptada es que los factores constantes en el tamiz de Atkin son mucho más altos y, por lo tanto, solo puede pagar por encima de los números de 32 bits ( N > 2 32 ) , si es así .

  

Editar: lo único que todavía no entiendo es a qué se refieren las variables x e y en el pseudocódigo. ¿Podría alguien arrojarme algo de luz sobre esto?

Para una explicación básica del uso común de las variables 'x' e 'y' en el seudocódigo de la página de Wikipedia (u otras implementaciones mejores del Tamiz de Atkin), puede encontrar mi respuesta a otra pregunta relacionada útil.

Aquí hay una implementación en C ++ de Tamiz de Atkins que podría ayudarlo ...

#include <iostream>
#include <cmath>
#include <fstream>
#include<stdio.h>
#include<conio.h>
using namespace std;

#define  limit  1000000

int root = (int)ceil(sqrt(limit));
bool sieve[limit];
int primes[(limit/2)+1];

int main (int argc, char* argv[])
{
   //Create the various different variables required
   FILE *fp=fopen("primes.txt","w");
   int insert = 2;
   primes[0] = 2;
   primes[1] = 3;
   for (int z = 0; z < limit; z++) sieve[z] = false; //Not all compilers have false as the       default boolean value
   for (int x = 1; x <= root; x++)
   {
        for (int y = 1; y <= root; y++)
        {
             //Main part of Sieve of Atkin
             int n = (4*x*x)+(y*y);
             if (n <= limit && (n % 12 == 1 || n % 12 == 5)) sieve[n] ^= true;
             n = (3*x*x)+(y*y);
             if (n <= limit && n % 12 == 7) sieve[n] ^= true;
             n = (3*x*x)-(y*y);
             if (x > y && n <= limit && n % 12 == 11) sieve[n] ^= true;
        }
   }
        //Mark all multiples of squares as non-prime
   for (int r = 5; r <= root; r++) if (sieve[r]) for (int i = r*r; i < limit; i += r*r) sieve[i] = false;
   //Add into prime array
   for (int a = 5; a < limit; a++)
   {
            if (sieve[a])
            {
                  primes[insert] = a;
                  insert++;
            }
   }
   //The following code just writes the array to a file
   for(int i=0;i<1000;i++){
             fprintf(fp,"%d",primes[i]);
             fprintf(fp,", ");
   }

   return 0;
 }
// Title : Seive of Atkin ( Prime number Generator) 

#include <iostream>
#include <cmath>
#include <vector>

using namespace std;

int main()
{
    ios_base::sync_with_stdio(false);
    long long int n;
    cout<<"Enter the value of n : ";
    cin>>n;
    vector<bool> is_prime(n+1);
    for(long long int i = 5; i <= n; i++)
    {
        is_prime[i] = false;
    }
    long long int lim = ceil(sqrt(n));

    for(long long int x = 1; x <= lim; x++)
    {
        for(long long int y = 1; y <= lim; y++)
        {
            long long int num = (4*x*x+y*y);
            if(num <= n && (num % 12 == 1 || num%12 == 5))
            {
                is_prime[num] = true;
            }

            num = (3*x*x + y*y);
            if(num <= n && (num % 12 == 7))
            {
                is_prime[num] = true;
            }

            if(x > y)
            {
                num = (3*x*x - y*y);
                if(num <= n && (num % 12 == 11))
                {
                    is_prime[num] = true;
                }
            }
        }
    }
    // Eliminating the composite by seiveing
    for(long long int i = 5; i <= lim; i++)
    {
        if(is_prime[i])
            for(long long int j = i*i; j <= n; j += i)
            {
                is_prime[j] = false;
            }
    }
    // Here we will start printing of prime number
   if(n > 2)
   {
       cout<<"2\t"<<"3\t";
   }
   for(long long int i = 5; i <= n; i++)
   {
         if(is_prime[i])
         {
             cout<<i<<"\t";
         }
    }
    cout<<"\n";
    return 0;
}
Licenciado bajo: CC-BY-SA con atribución
No afiliado a StackOverflow
scroll top