O algoritmo mais conhecido para encontrar todos os números primos não maiores que um dado é a peneira de Eratóstenes . Funciona muito bem para números de até bilhões, talvez até dezenas de bilhões, se bem escritos. No entanto, todo mundo que gosta de se divertir com números primos sabe que você sempre quer ter o maior número possível em mãos. Uma vez, para resolver um problema em um hacker, eu precisava de um banco de dados na memória de números primos de até cem bilhões. Com a otimização máxima da memória, se você representar números ímpares em uma peneira de Eratóstenes com uma matriz de bits, seu tamanho será de cerca de 6 gigabytes, o que não coube na memória do meu laptop. Há uma modificação do algoritmo que é muito menos exigente na memória (dividindo o intervalo original de números em várias partes e processando uma parte de cada vez) -peneira segmentada de Eratóstenes , mas é mais difícil de implementar, e todo o resultado não se encaixará na memória. Abaixo, chamo a atenção um algoritmo que é quase tão simples quanto a peneira de Eratóstenes, mas que oferece otimização dupla pela memória (ou seja, um banco de dados de primos de até cem bilhões ocupará cerca de 3 gigabytes, que já devem caber na memória de um laptop padrão).Primeiro, vamos relembrar como a peneira de Eratóstenes funciona: em uma matriz de números de 1 a N, números divisíveis por dois são cruzados, depois números divisíveis por 3 e assim por diante. Os números restantes após o cruzamento serão simples:2 3 . 5 . 7 . 9 . 11 . 13 . 15 . 17 . 19 . 21 . 23 . 25 . 27 . 29 ...
2 3 . 5 . 7 . . . 11 . 13 . . . 17 . 19 . . . 23 . 25 . . . 29 ...
2 3 . 5 . 7 . . . 11 . 13 . . . 17 . 19 . . . 23 . . . . . 29 ...
Código c #class EratospheneSieve
{
private bool[] Data;
public EratospheneSieve(int length)
{
Data = new bool[length];
for (int i = 2; i < Data.Length; i++) Data[i] = true;
for (int i = 2; i * i < Length; i++)
{
if (!Data[i]) continue;
for (int d = i * i; d < Length; d += i) Data[d] = false;
}
}
public int Length => Data.Length;
public bool IsPrime(int n) => Data[n];
}
No meu laptop, ele conta até 15 bilhões de segundos e consome um gigabyte de memória.O algoritmo que usaremos é chamado de fatoração de roda . Para entender sua essência, escrevemos números naturais com um tablet de 30 peças seguidas: 0 1 2 3 4 5 ... 26 27 28 29
30 31 32 33 34 35 ... 56 57 58 59
60 61 62 63 64 65 ... 86 87 88 89
...
E cruze todos os números divisíveis por 2, 3 e 5:. 1 . . . . . 7 . . . 11 . 13 . . . 17 . 19 . . . 23 . . . . . 29
. 31 . . . . . 37 . . . 41 . 43 . . . 47 . 49 . . . 53 . . . . . 59
. 61 . . . . . 67 . . . 71 . 73 . . . 77 . 79 . . . 83 . . . . . 89
...
Pode-se ver que em cada linha os números são riscados nas mesmas posições. Como cruzamos apenas números compostos (com a exceção, de fato, dos números 2, 3, 5), os números primos só podem estar em posições não marcadas.Como existem oito números em cada linha, isso nos dá a idéia de representar cada trinta números com um byte, sendo que cada um dos oito bits indicará a simplicidade do número na posição correspondente. Essa coincidência não é apenas esteticamente bonita, mas também simplifica bastante a implementação do algoritmo!
Implementação técnica
Primeiro, nosso objetivo é atingir cem bilhões, para que os números de quatro bytes não possam mais ser dispensados. Portanto, nosso algoritmo será mais ou menos assim:class Wheel235
{
private byte[] Data;
public Wheel235(long length)
{
...
}
public bool IsPrime(long n)
{
...
}
...
}
Para facilitar a implementação, arredondaremos o comprimento até o número mais próximo divisível por 30. Em seguida,class Wheel235
{
private byte[] Data;
public Wheel235(long length)
{
length = (length + 29) / 30 * 30;
Data = new byte[length/30];
...
}
public long Length => Data.Length * 30L;
...
}
Em seguida, precisamos de duas matrizes: uma traduz o índice 0..29 no número do bit, a segunda, pelo contrário, o número do bit na posição trinta. Usando essas matrizes e um pouco de aritmética de bits, construímos uma correspondência individual entre números naturais e bits na matriz de dados. Para simplificar, usamos apenas dois métodos: IsPrime (n), que permite descobrir se o número é primo, e ClearPrime (n), que define o bit correspondente como 0, marcando o número n como composto:class Wheel235
{
private static int[] BIT_TO_INDEX = new int[] { 1, 7, 11, 13, 17, 19, 23, 29 };
private static int[] INDEX_TO_BIT = new int[] {
-1, 0,
-1, -1, -1, -1, -1, 1,
-1, -1, -1, 2,
-1, 3,
-1, -1, -1, 4,
-1, 5,
-1, -1, -1, 6,
-1, -1, -1, -1, -1, 7,
};
private byte[] Data;
...
public bool IsPrime(long n)
{
if (n <= 5) return n == 2 || n == 3 || n == 5;
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return false;
return (Data[n / 30] & (1 << bit)) != 0;
}
private void ClearPrime(long n)
{
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return;
Data[n / 30] &= (byte)~(1 << bit);
}
}
Observe que se o número for divisível por 2, 3 ou 5, o número do bit será -1, o que significa que o número é obviamente composto. Bem, e claro, estamos processando um caso especial n <= 5.Agora, de fato, o próprio algoritmo, que quase literalmente repete a peneira de Eratóstenes: public Wheel235(long length)
{
length = (length + 29) / 30 * 30;
Data = new byte[length / 30];
for (long i = 0; i < Data.Length; i++) Data[i] = byte.MaxValue;
for (long i = 7; i * i < Length; i++)
{
if (!IsPrime(i)) continue;
for (long d = i * i; d < Length; d += i) ClearPrime(d);
}
}
Assim, toda a classe resultante:parece queclass Wheel235
{
private static int[] BIT_TO_INDEX = new int[] { 1, 7, 11, 13, 17, 19, 23, 29 };
private static int[] INDEX_TO_BIT = new int[] {
-1, 0,
-1, -1, -1, -1, -1, 1,
-1, -1, -1, 2,
-1, 3,
-1, -1, -1, 4,
-1, 5,
-1, -1, -1, 6,
-1, -1, -1, -1, -1, 7,
};
private byte[] Data;
public Wheel235(long length)
{
length = (length + 29) / 30 * 30;
Data = new byte[length / 30];
for (long i = 0; i < Data.Length; i++) Data[i] = byte.MaxValue;
for (long i = 7; i * i < Length; i++)
{
if (!IsPrime(i)) continue;
for (long d = i * i; d < Length; d += i) ClearPrime(d);
}
}
public long Length => Data.Length * 30L;
public bool IsPrime(long n)
{
if (n >= Length) throw new ArgumentException("Number too big");
if (n <= 5) return n == 2 || n == 3 || n == 5;
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return false;
return (Data[n / 30] & (1 << bit)) != 0;
}
private void ClearPrime(long n)
{
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return;
Data[n / 30] &= (byte)~(1 << bit);
}
}
Há apenas um pequeno problema: esse código não funciona para números maiores que cerca de 65 bilhões! Quando você tenta iniciá-lo com esses números, o programa trava com um erro:Unhandled Exception: System.OverflowException: Array dimensions exceeded supported range.
O problema é que em matrizes c # não pode ter mais de 2 ^ 31 elementos (aparentemente porque a implementação interna usa um índice de matriz de quatro bytes). Uma opção é criar, por exemplo, uma matriz longa [] em vez de uma matriz de bytes [], mas isso complicará um pouco a aritmética dos bits. Para simplificar, seguiremos o outro caminho, criando uma classe especial que simula a matriz desejada, mantendo duas matrizes curtas dentro. Ao mesmo tempo, daremos a ele a oportunidade de nos salvar em disco, para que você possa calcular os primos uma vez e, em seguida, use o banco de dados novamente: class LongArray
{
const long MAX_CHUNK_LENGTH = 2_000_000_000L;
private byte[] firstBuffer;
private byte[] secondBuffer;
public LongArray(long length)
{
firstBuffer = new byte[Math.Min(MAX_CHUNK_LENGTH, length)];
if (length > MAX_CHUNK_LENGTH) secondBuffer = new byte[length - MAX_CHUNK_LENGTH];
}
public LongArray(string file)
{
...
}
public long Length => firstBuffer.LongLength + (secondBuffer == null ? 0 : secondBuffer.LongLength);
public byte this[long index]
{
get => index < MAX_CHUNK_LENGTH ? firstBuffer[index] : secondBuffer[index - MAX_CHUNK_LENGTH];
set
{
if (index < MAX_CHUNK_LENGTH) firstBuffer[index] = value;
else secondBuffer[index - MAX_CHUNK_LENGTH] = value;
}
}
public void Save(string file)
{
...
}
}
Por fim, combine todas as etapas emversão final do algoritmoclass Wheel235
{
class LongArray
{
const long MAX_CHUNK_LENGTH = 2_000_000_000L;
private byte[] firstBuffer;
private byte[] secondBuffer;
public LongArray(long length)
{
firstBuffer = new byte[Math.Min(MAX_CHUNK_LENGTH, length)];
if (length > MAX_CHUNK_LENGTH) secondBuffer = new byte[length - MAX_CHUNK_LENGTH];
}
public LongArray(string file)
{
using(FileStream stream = File.OpenRead(file))
{
long length = stream.Length;
firstBuffer = new byte[Math.Min(MAX_CHUNK_LENGTH, length)];
if (length > MAX_CHUNK_LENGTH) secondBuffer = new byte[length - MAX_CHUNK_LENGTH];
stream.Read(firstBuffer, 0, firstBuffer.Length);
if(secondBuffer != null) stream.Read(secondBuffer, 0, secondBuffer.Length);
}
}
public long Length => firstBuffer.LongLength + (secondBuffer == null ? 0 : secondBuffer.LongLength);
public byte this[long index]
{
get => index < MAX_CHUNK_LENGTH ? firstBuffer[index] : secondBuffer[index - MAX_CHUNK_LENGTH];
set
{
if (index < MAX_CHUNK_LENGTH) firstBuffer[index] = value;
else secondBuffer[index - MAX_CHUNK_LENGTH] = value;
}
}
public void Save(string file)
{
using(FileStream stream = File.OpenWrite(file))
{
stream.Write(firstBuffer, 0, firstBuffer.Length);
if (secondBuffer != null) stream.Write(secondBuffer, 0, secondBuffer.Length);
}
}
}
private static int[] BIT_TO_INDEX = new int[] { 1, 7, 11, 13, 17, 19, 23, 29 };
private static int[] INDEX_TO_BIT = new int[] {
-1, 0,
-1, -1, -1, -1, -1, 1,
-1, -1, -1, 2,
-1, 3,
-1, -1, -1, 4,
-1, 5,
-1, -1, -1, 6,
-1, -1, -1, -1, -1, 7,
};
private LongArray Data;
public Wheel235(long length)
{
length = (length + 29) / 30 * 30;
Data = new LongArray(length / 30);
for (long i = 0; i < Data.Length; i++) Data[i] = byte.MaxValue;
for (long i = 7; i * i < Length; i++)
{
if (!IsPrime(i)) continue;
for (long d = i * i; d < Length; d += i) ClearPrime(d);
}
}
public Wheel235(string file)
{
Data = new LongArray(file);
}
public void Save(string file) => Data.Save(file);
public long Length => Data.Length * 30L;
public bool IsPrime(long n)
{
if (n >= Length) throw new ArgumentException("Number too big");
if (n <= 5) return n == 2 || n == 3 || n == 5;
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return false;
return (Data[n / 30] & (1 << bit)) != 0;
}
private void ClearPrime(long n)
{
int bit = INDEX_TO_BIT[n % 30];
if (bit < 0) return;
Data[n / 30] &= (byte)~(1 << bit);
}
}
Como resultado, obtivemos uma classe que pode calcular números primos não superiores a cem bilhões e salvar o resultado em disco (eu tinha esse código no meu laptop por 50 minutos; o tamanho do arquivo criado era de cerca de 3 gigabytes, conforme o esperado): long N = 100_000_000_000L;
Wheel235 wheel = new Wheel235(N);
wheel.Save("./primes.dat");
E depois leia a partir do disco e use: DateTime start = DateTime.UtcNow;
Wheel235 loaded = new Wheel235("./primes.dat");
DateTime end = DateTime.UtcNow;
Console.WriteLine($"Database of prime numbers up to {loaded.Length:N0} loaded from file to memory in {end - start}");
long number = 98_000_000_093L;
Console.WriteLine($"{number:N0} is {(loaded.IsPrime(number) ? "prime" : "not prime")}!");