Реализации алгоритмов/Решето Аткина: различия между версиями

Содержимое удалено Содержимое добавлено
Добавлены ссылки на Решето Сундарама и Эратосфена
Оптимизация кода
Строка 19:
* Флаги «простое\составное» хранятся упакованными в массив битов (<code>vector<bool></code> для С++, <code>BitSet</code> для Java и т. п.). Значение ''n''-го бита ''true'' соответствует простому ''n'', ''false'' — составному.
* Умножение на степень двойки производится операцией побитового сдвига (''<<''). Например <code>a << 2</code> вместо <code>4 * a</code>.
* Для уменьшения риска переполнения вычисления ведутся в типе длинное целое (<code>unsigned long long</code> в C++, <code>long</code> в Java); на выходе значения выдаются в заданном целочисленном типе.
* Квадраты последовательных натуральных чисел вычисляются по формуле ''(x + 1)² = x² + 2x + 1'', т. е. прибавлением последовательных нечётных чисел. Другими словами, если ''squares'' — последовательность квадратов натуральных чисел ''{0² = 0, 1² = 1, 2² = 4…}'', а ''odds'' — последовательность нечётных чисел ''{1, 3, 5…}'', то ''square<sub>i+1</sub> = square<sub>i</sub> + odds<sub>i</sub>'' (для любого ''i ≥ 0'').
* ''n'' делится с остатком на 12 (а не на 60). Полученные значения 1 и 5 соответствуют случаю (1), 7 — (2), 11 — (3), прочие — составному ''n''.
Строка 26 ⟶ 27 :
 
<source lang="cpp">
// atkinatkin_sieve.hpp
#include <vector>
 
 
namespace AtkinAtkinSieve {
template < typename T, typename Sequence >
namespace {
void sievegetPrimesUpTo (Tunsigned const limit, Sequence &primes) {
typedef std::vector<bool> _Sieve;
}
 
template < typename T, typename Sequence >
void sieve (T const limit, Sequence &primes) {
// Очистка результирующей последовательности
primes.clear();
// Числа 2 и 3 — заведомо простые
if (limit > 22u)
primes.push_back(22u);
else
return;
if (limit > 33u)
primes.push_back(33u);
else
return;
// Предварительное просеивание
Tunsigned long long n;
_Sievestd::vector<bool> sieve(limit);
for (Tunsigned long long x2 = 11ull, dx2 = 33ull; x2 < limit; x2 += dx2, dx2 += 22ull)
for (Tunsigned long long y2 = 11ull, dy2 = 33ull; y2 < limit; y2 += dy2, dy2 += 22ull) {
// n = 4x² + y²
n = (x2 << 22ull) + y2;
if (n <= limit && (n % 1212ull == 11ull || n % 1212ull == 55ull))
sieve[n].flip();
// n = 3x² + y²
n -= x2;
if (n <= limit && n % 1212ull == 77ull)
sieve[n].flip();
// n = 3x² - y² (при x > y)
if (x2 > y2) {
n -= y2 << 11ull;
if (n <= limit && n % 1212ull == 1111ull)
sieve[n].flip();
}
}
// Все числа, кратные квадратам, помечаются как составные
for (Tunsigned long long r = 55ull, r2 = r * r, dr2 = (r << 11ull) + 11ull; r2 < limit; ++r, r2 += dr2, dr2 += 22ull)
if (sieve[r])
for (n = r2; n < limit; n += r2)
sieve[n] = false;
// Занесение чисел в результирующую последовательность
for (n = 55ull; n < limit; ++n)
if (sieve[n])
primes.push_back(n);
}
}</source>
}
 
</source>
Пример использования:
<source lang="cpp">
Строка 86 ⟶ 83 :
#include <vector>
 
#include "atkinatkin_sieve.hpp"
 
 
Строка 93 ⟶ 90 :
int main () {
Sequence primes;
AtkinAtkinSieve::sievegetPrimesUpTo(100u, primes); // Вычисление простых чисел ≤ 100
// Вывод последовательности
for (Sequence::const_iterator iterator = primes.begin(); iterator != primes.end(); ++iterator)
Строка 104 ⟶ 101 :
Простые числа вычисляются для типа <code>int</code>. Результирующая последовательность должна иметь тип, реализующий интерфейс <code>List<Integer></code> (например, <code>ArrayList<Integer></code>).
<source lang="java">
// AtkinAtkinSieve.java
import java.util.BitSet;
import java.util.List;
 
 
public class AtkinAtkinSieve {
public static void sievegetPrimesUpTo (int limit, List<Integer> primes) {
// Очистка результирующей последовательности
primes.clear();
Строка 123 ⟶ 120 :
return;
// Предварительное просеивание
intlong n;
BitSet sieve = new BitSet();
for (intlong x2 = 11L, dx2 = 33L; x2 < limit; x2 += dx2, dx2 += 22L)
for (intlong y2 = 11L, dy2 = 33L; y2 < limit; y2 += dy2, dy2 += 22L) {
// n = 4x² + y²
n = (x2 << 22L) + y2;
if (n <= limit && (n % 1212L == 11L || n % 1212L == 55L))
sieve.flip((int)n);
// n = 3x² + y²
n -= x2;
if (n <= limit && n % 1212L == 77L)
sieve.flip((int)n);
// n = 3x² - y² (при x > y)
if (x2 > y2) {
n -= y2 << 11L;
if (n <= limit && n % 1212L == 1111L)
sieve.flip((int)n);
}
}
// Все числа, кратные квадратам, помечаются как составные
for (intlong r = 55L, r2 = r * r, dr2 = (r << 11L) + 11L; r2 < limit; ++r, r2 += dr2, dr2 += 22L)
if (sieve.get(r))
for (n = r2; n < limit; n += r2)
sieve.set((int)n, false);
// Занесение чисел в результирующую последовательность
for (nint i = 55L; ni < limit; ++ni)
if (sieve.get(n)i)
primes.add(ni);
}
}
Строка 163 ⟶ 160 :
public static void main (String []args){
ArrayList<Integer> primes = new ArrayList<>();
AtkinAtkinSieve.sievegetPrimesUpTo(100, primes); // Вычисление простых чисел ≤ 100
// Вывод последовательности
for (int prime : primes)