quarta-feira, 12 de janeiro de 2022

O crivo de Eratóstenes

Segundo  Coutinho (2005): "o crivo de Eratóstenes é o mais antigo dos métodos para achar primos, e não envolve nenhuma fórmula explícita." Crivo significa peneira e o método consiste em "passar" os números por uma peneira, sobrando apenas os números primos. Dado um número inteiro n positivo, o crivo de Eratóstenes encontra todos os números primos menores que n. Vamos entender como esse método funciona. (lembrando que qualquer número par maior que 2 é composto).

Listamos os números ímpares de 3 a n. Começamos então a riscar os números múltiplos de 3, ou seja, a partir do número 3 riscamos números de 3 em 3. Fazemos o mesmo passo para o número 5, ou seja, riscamos números de 5 em 5 a partir de 5. E assim sucessivamente, até o último número da nossa lista. Os números não riscados ao final do método são números primos.

Vamos exemplificar. Suponha que listamos os números ímpares de 3 até 29.

3 5 7 9 11 13 15 17 19 21 23 25 27 29

A primeira etapa é pegar o primeiro número da lista, que é o 3, e riscar os números múltiplos de 3:

3 5 7 9 11 13 15 17 19 21 23 25 27 29

Pegamos o próximo número não riscado da lista, que é o 5, e riscamos os números múltiplos de 5:

3 5 7 9 11 13 15 17 19 21 23 25 27 29

O processo continua até acabarem os números não riscados. O resultado final do nosso exemplo é:

3 5 7 9 11 13 15 17 19 21 23 25 27 29

O que significa que os números 3, 5, 7, 11, 13, 17, 19, 23 e 29 são primos. Perceba que não realizamos contas nesse método, apenas riscamos os números múltiplos de um certo número em uma determinada etapa do método.

Agora, vamos ao que interessa: implementar o método!

Vamos começar definindo o vetor de números ímpares:


n = 29
v_num = seq(3,n,by=2)

> v_num
 [1]  3  5  7  9 11 13 15 17 19 21 23 25 27 29

Para sabermos quais números foram riscados e quais não foram, vamos definir um vetor do mesmo tamanho que v_num, cujas entradas sejam todas iguais a 1, ou seja, inicialmente não temos nenhum número riscado. O valor zero vai indicar que o número foi riscado.


tam_v = length(v_num)
v_i = rep(1,tam_v)
rbind(v_num,v_i)

> rbind(v_num,v_i)
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14]
v_num    3    5    7    9   11   13   15   17   19    21    23    25    27    29
v_i      1    1    1    1    1    1    1    1    1     1     1     1     1     1

Desse modo, v_num[1] armazena o número 3 e v_i[1] armazena se o número 3 foi riscado (0) ou não (1).

Agora entramos na implementação do método:


Para cada elemento do vetor v_num
	pega o número primo atual
  	seja aux = indice atual + número primo atual
    	enquanto aux for menor que o tamanho do vetor v_num
    		altera para 0 o vetor v_i no indice aux
        	incrementa aux do valor do número primo atual
        fim enquanto
fim do para

O pseudo código acima pode ser escrito da seguinte forma no R:


for(i in 1:tam_v){
  p_atual = v_num[i];
  aux = i + p_atual;
  while(aux < tam_v){
    v_i[aux] = 0;
    aux = aux + p_atual;
  }
}

Ou seja, para cada elemento do vetor v_num, armazenamos o número primo atual em p_atual, armazenamos o próximo múltiplo desse número na variável aux e, enquanto não percorremos todo o vetor de números ímpares, vamos alterando para 0 o valor de v_i no índice correspondente aos múltiplos do primo atual. Ao final do processo, temos:


> rbind(v_num,v_i)
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [,11] [,12] [,13] [,14]
v_num    3    5    7    9   11   13   15   17   19    21    23    25    27    29
v_i      1    1    1    0    1    1    0    1    1     0     1     0     0     1

Note que o resultado acima é exatamente o resultado do crivo de Eratóstenes: os números primos não estão riscados (1) e os números compostos estão riscados (0). Para selecionarmos apenas os número primos, fazemos:


> v_num[which(v_i==1)]
[1]  3  5  7 11 13 17 19 23 29

Como o código está solto, vamos colocá-lo dentro de uma função:


crivo_erat = function(n){
  v_num = seq(3,n,by=2)
  tam_v = length(v_num)
  v_i = rep(1,tam_v)
  
  for(i in 1:tam_v){
    p_atual = v_num[i];
    aux = i + p_atual;
    while(aux < tam_v){
      v_i[aux] = 0;
      aux = aux + p_atual;
    }
    
  }
  primos = v_num[which(v_i==1)];
  
  return(primos)
}

Desse modo, nossa função recebe um argumento: o número n escolhido, e retorna o resultado: um vetor contendo os números primos. Vamos testar:


> crivo_erat(100)
 [1]  3  5  7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97 99
 
> crivo_erat(1000)
  [1]   3   5   7  11  13  17  19  23  29  31  37  41  43  47  53  59  61  67  71  73  79  83
 [23]  89  97 101 103 107 109 113 127 131 137 139 149 151 157 163 167 173 179 181 191 193 197
 [45] 199 211 223 227 229 233 239 241 251 257 263 269 271 277 281 283 293 307 311 313 317 331
 [67] 337 347 349 353 359 367 373 379 383 389 397 401 409 419 421 431 433 439 443 449 457 461
 [89] 463 467 479 487 491 499 503 509 521 523 541 547 557 563 569 571 577 587 593 599 601 607
[111] 613 617 619 631 641 643 647 653 659 661 673 677 683 691 701 709 719 727 733 739 743 751
[133] 757 761 769 773 787 797 809 811 821 823 827 829 839 853 857 859 863 877 881 883 887 907
[155] 911 919 929 937 941 947 953 967 971 977 983 991 997 999

Agora que temos o crivo de Eratóstenes implementado, podemos apontar dois pontos ineficientes da nossa implementação:

1) Alguns números são riscados da lista mais de uma vez.

2) Não é necessário percorrer todo o vetor de números ímpares para obtermos o resultado.

Vamos discutir esses dois pontos com o exemplo de uma lista de números de 3 a 49:

3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49

Primeiro ponto: na etapa do número primo 5, o número 15 já foi riscado na etapa do número primo 3. Note que 15 = 5*3, ou seja, o próximo número não riscado múltiplo de 5 do nosso vetor de números ímpares seria o número 5*5 = 25. Na etapa do número primo 7, o número 21 = 7*3 foi riscado por ser múltiplo de 3 e o número 35 = 7*5 foi riscado por ser múltiplo de 5. Ou seja, o próximo número não riscado múltiplo de 7 é 7*7 = 49. Dessa forma, dado um número primo atual p, podemos começar a riscar números a partir de p^2, evitando cortes desnecessários e tornando o algoritmo mais eficiente.

Segundo ponto: podemos encerrar o algoritmo quando o número primo atual for maior que a raiz quadrada do maior número da lista de ímpares n. Por quê? Lembrem das aulas passadas, quando otimizamos nosso algoritmo para encontrar números primos: um número composto terá um fator menor ou igual a sua raiz quadrada. Logo, se quisermos saber se 49 é número primo, basta verificarmos divisão até o número 7. Mas e os demais números do crivo? Esses números serão menores que o maior número da lista e, portanto, a verificação é ainda menor. Por exemplo, paramos a verificação do número 47 quando o primo atual é maior que 6,8. Isso significa que, depois do número 7, não cortaremos nada da lista de números ímpares.
Não te convenci? Vamos ao exemplo.
Após cortamos múltiplos de 3:

3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49

Após cortamos múltiplos de 5:

3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49

Após cortamos múltiplos de 7:

3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49

Agora, na etapa do número 11, o número 33 já foi cortado e não há mais múltiplos de 11 na nossa lista.

Na etapa do número 13, o número 39 já foi cortado e não há mais múltiplos de 13 na nossa lista.

Perceba que, se todos os múltiplos do número atual menores que seu quadrado já foram cortados, não há mais o que fazer nesse exemplo. No caso do número 11, começaríamos a cortar números a partir de 121, que não está na nossa lista. Dessa forma, podemos encerrar o algoritmo na etapa do número 7, que é a raiz quadrada de maior número da nossa lista.

Vamos implementar esses dois pontos no nosso algoritmo:


crivo_erat2 = function(n){
  v_num = seq(3,n,by=2)
  tam_v = length(v_num)
  v_i = rep(1,tam_v)
  ult_num = v_num[tam_v]
  
  for(i in 1:tam_v){
    if(v_num[i]>sqrt(ult_num)) break;
    p_atual = v_num[i];
    aux = (p_atual^2 - 1)/2;
    while(aux < tam_v){
      v_i[aux] = 0;
      aux = aux + p_atual;
    }
  }
  primos = v_num[which(v_i==1)];
  
  return(primos)
}

Vejamos as alterações que fizemos. O primeiro ponto está na seguinte linha de código


    aux = (p_atual^2 - 1)/2;

A linha de código acima diz o seguinte: o primeiro índice que pegaremos será o correspondente ao quadrado do número primo atual. Por exemplo, se estamos na etapa no número 7, começamos os cortes a partir do número 49, que tem como índice (49-1)/2 = 24 na nossa lista de ímpares.

O segundo ponto está implementado na seguinte linha de código


    if(v_num[i]>sqrt(ult_num)) break;

Simplesmente paramos o for quando o número primo atual é maior que a raiz quadrada do maior número da nossa lista.

Beleza, implementamos as melhorias, ou otimizações. Mas você deve estar se perguntando: vale a pena todo esse trabalho? É claro que a implementação não deu trabalho, são linhas de código simples, mas pensá-las é o que realmente dá trabalho e nem sempre é simples. Será que vale o esforço? Vamos analisar o custo de tempo das duas funções utilizando a função benchmark do pacote rbenchmark. Vamos aproveitar e comparar o tempo de execução com a função calcula_primos2() que implementamos nas aulas passadas.


> numero = 200;
> benchmark(crivo_erat(numero),crivo_erat2(numero),
+           calcula_primos2(numero),
+           replications = 1000,
+           columns = c("test", "replications", "elapsed","relative"))
                     test replications elapsed relative
3 calcula_primos2(numero)         1000    0.24      4.8
1      crivo_erat(numero)         1000    0.06      1.2
2     crivo_erat2(numero)         1000    0.05      1.0

> numero = 2000;
> benchmark(crivo_erat(numero),crivo_erat2(numero),
+           calcula_primos2(numero),
+           replications = 1000,
+           columns = c("test", "replications", "elapsed","relative"))
                     test replications elapsed relative
3 calcula_primos2(numero)         1000    3.06   16.105
1      crivo_erat(numero)         1000    0.56    2.947
2     crivo_erat2(numero)         1000    0.19    1.000

Percebam que sim! valeu a pena implementar as melhorias do crivo no nosso algoritmo. Vejam que mesmo para n pequeno, a função crivo_erat2 foi a de menor tempo. Para n = 2000, vemos uma grande diferença de tempo entre as duas funções do crivo, em que a função não otimizada leva quase 3 vezes o tempo de execução da função otimizada.

Agora, em comparação com a nossa função calcula_primos2 das aulas passadas (que já está otimizada), percebam que o crivo de Eratóstenes é bem mais eficiente. O que isso significa? Que apenas otimização de uma função não garante que ela seja a melhor, mas sim que ela é melhor que sua versão não otimizada. Podem existir métodos mais eficientes do qual você esteja usando, que valem mais a pena serem implementados do que simplesmente otimizar a função do método atual. Lembrem-se que, na prática, apesar de as funções (e métodos) retornarem o mesmo resultado, o melhor é escolher o mais eficiente!

Até a próxima aula!

Referência citada:

Coutinho, S.C. Números Inteiros e Criptografia RSA. 2005. Série de Computação e Matemática. IMPA.

Nenhum comentário:

Postar um comentário