Эратосфены шигшүүр¶

Эратосфены шигшүүр бол $[1;n]$ хэрчим дэх бүх анхны тоог $O(n \log \log n)$ үйлдлээр олох алгоритм юм.

Алгоритм нь маш энгийн: эхлээд бид 2-оос $n$ хүртэлх бүх тоог бичнэ. 2-ын бүх жинхэнэ үржвэрийг (2 бол хамгийн бага анхны тоо тул) нийлмэл гэж тэмдэглэнэ. $x$ тооны жинхэнэ үржвэр гэдэг нь $x$-ээс их бөгөөд $x$-д хуваагддаг тоо юм. Дараа нь нийлмэл гэж тэмдэглэгдээгүй дараагийн тоог олно, энэ тохиолдолд 3 болно. Энэ нь 3 анхны тоо гэсэн үг бөгөөд бид 3-ын бүх жинхэнэ үржвэрийг нийлмэл гэж тэмдэглэнэ. Дараагийн тэмдэглэгдээгүй тоо нь 5 бөгөөд энэ бол дараагийн анхны тоо, түүний бүх жинхэнэ үржвэрийг мөн тэмдэглэнэ. Мөрөн дэх бүх тоог боловсруулж дуустал бид энэ үйл ажиллагааг үргэлжлүүлнэ.

Дараах зурган дээр $[1; 16]$ мужийн бүх анхны тоог олох алгоритмын дүрслэлийг харж болно. Бид тоонуудыг нэлээд олон удаа давхардуулан нийлмэл гэж тэмдэглэдгийг харж болно.

Үүний ард байгаа санаа нь: Хэрэв түүнээс бага ямар ч анхны тоо хуваахгүй бол тухайн тоо анхны тоо болно. Бид анхны тоонуудыг эрэмбээр давтаж байгаа тул дор хаяж нэг анхны тоонд хуваагддаг бүх тоог аль хэдийн хуваагдана гэж тэмдэглэсэн байна. Тиймээс бид нүд рүү хүрэхэд тэмдэглэгдээгүй байвал энэ нь өөрөөсөө бага ямар ч анхны тоонд хуваагдахгүй тул анхны тоо байх ёстой.

Implementation¶

int n;
vector<bool> is_prime(n+1, true);
is_prime[0] = is_prime[1] = false;
for (int i = 2; i <= n; i++) {
    if (is_prime[i] && (long long)i * i <= n) {
        for (int j = i * i; j <= n; j += i)
            is_prime[j] = false;
    }
}

This code first marks all numbers except zero and one as potential prime numbers, then it begins the process of sifting composite numbers. For this it iterates over all numbers from $2$ to $n$. If the current number $i$ is a prime number, it marks all numbers that are multiples of $i$ as composite numbers, starting from $i^2$. This is already an optimization over naive way of implementing it, and is allowed as all smaller numbers that are multiples of $i$ necessary also have a prime factor which is less than $i$, so all of them were already sifted earlier. Since $i^2$ can easily overflow the type int, the additional verification is done using type long long before the second nested loop.

Using such implementation the algorithm consumes $O(n)$ of the memory (obviously) and performs $O(n \log \log n)$ (see next section).

Асимптот шинжилгээ¶

Анхны тоонуудын тархалтын талаар юу ч мэдэлгүйгээр $O(n \log n)$ ажиллах хугацааг батлахад амархан — is_prime шалгалтыг үл тоомсорловол дотоод давталт $i = 2, 3, 4, \dots$-ийн хувьд (хамгийн ихдээ) $n/i$ удаа ажиллах бөгөөд ингэснээр дотоод давталтын нийт үйлдлийн тоо $n(1/2 + 1/3 + 1/4 + \cdots)$ гэсэн гармоник нийлбэр болж, $O(n \log n)$-ээр хязгаарлагдана.

Алгоритмын ажиллах хугацаа $O(n \log \log n)$ болохыг баталъя. Алгоритм нь дотоод давталтад $n$-ээс бага буюу тэнцүү анхны тоо $p$ бүрийн хувьд $\frac{n}{p}$ үйлдэл хийнэ. Тиймээс бид дараах илэрхийллийг үнэлэх хэрэгтэй:

$$\sum_{\substack{p \le n, \\\ p \text{ prime}}} \frac n p = n \cdot \sum_{\substack{p \le n, \\\ p \text{ prime}}} \frac 1 p.$$

Мэдэгдэж буй хоёр баримтыг эргэн санацгаая.

$n$-ээс бага буюу тэнцүү анхны тооны тоо ойролцоогоор $\frac n {\ln n}$ байна.
$k$ дугаар анхны тоо ойролцоогоор $k \ln k$-тэй тэнцүү (энэ нь өмнөх баримтаас гарна).

Ингэснээр бид нийлбэрийг дараах байдлаар бичиж болно:

$$\sum_{\substack{p \le n, \\\ p \text{ prime}}} \frac 1 p \approx \frac 1 2 + \sum_{k = 2}^{\frac n {\ln n}} \frac 1 {k \ln k}.$$

Энд бид эхний анхны тоо 2-ыг нийлбэрээс салгаж авсан, учир нь $k \ln k$ ойролцоололд $k = 1$ үед $0$ болж, тэгд хуваах алдаа үүсгэнэ.

Одоо энэ нийлбэрийг $k$-гаар $2$-оос $\frac n {\ln n}$ хүртэлх ижил функцийн интеграл ашиглан үнэлье (ийм ойролцоолол хийж болно, учир нь үнэндээ нийлбэр нь тэгш өнцөгтийн аргаар интегралыг ойролцоолсонтой холбоотой):

$$\sum_{k = 2}^{\frac n {\ln n}} \frac 1 {k \ln k} \approx \int_2^{\frac n {\ln n}} \frac 1 {k \ln k} dk.$$

Интеграл доорх функцийн эх функц нь $\ln \ln k$ юм. Орлуулга хийж, доод эрэмбийн гишүүдийг хасвал үр дүн нь:

$$\int_2^{\frac n {\ln n}} \frac 1 {k \ln k} dk = \ln \ln \frac n {\ln n} - \ln \ln 2 = \ln(\ln n - \ln \ln n) - \ln \ln 2 \approx \ln \ln n.$$

Одоо анхны нийлбэр рүүгээ буцаж очвол түүний ойролцоо үнэлгээг олно:

$$\sum_{\substack{p \le n, \\\ p\ is\ prime}} \frac n p \approx n \ln \ln n + o(n).$$

Илүү нарийн баталгааг (тогтмол үржүүлэгчийн нарийвчлалтай илүү оновчтой үнэлгээ өгдөг) Харди, Райт нарын бичсэн "An Introduction to the Theory of Numbers" номноос (х. 349) олж болно.

Эратосфены шигшүүрийн янз бүрийн оновчлол¶

Алгоритмын хамгийн том сул тал нь санах ойн дагуу олон удаа "алхаж", зөвхөн ганц элементтэй ажилладаг явдал юм. Энэ нь кэшийн хувьд тааламжгүй. Тиймээс $O(n \log \log n)$-д нуугдаж буй тогтмол харьцангуй том байдаг.

Түүнчлэн зарцуулагдах санах ой нь том $n$-ийн хувьд саад тотгор болдог.

Доор үзүүлсэн аргууд нь хийгдэх үйлдлийн тоог багасгах, мөн зарцуулагдах санах ойг мэдэгдэхүйц багасгах боломж олгодог.

Язгуур хүртэл шигших¶

$n$ хүртэлх бүх анхны тоог олохын тулд зөвхөн $n$-ийн язгуураас хэтрэхгүй анхны тоонуудаар шигших нь хангалттай нь ойлгомжтой.

int n;
vector<bool> is_prime(n+1, true);
is_prime[0] = is_prime[1] = false;
for (int i = 2; i * i <= n; i++) {
    if (is_prime[i]) {
        for (int j = i * i; j <= n; j += i)
            is_prime[j] = false;
    }
}

Ийм оновчлол нь complexity-д нөлөөлөхгүй (үнэндээ дээрх баталгааг давтвал $n \ln \ln \sqrt n + o(n)$ үнэлгээ гарах бөгөөд логарифмын шинж чанараар асимптотын хувьд ижил байна), гэвч үйлдлийн тоо мэдэгдэхүйц багасна.

Зөвхөн сондгой тоогоор шигших¶

Бүх тэгш тоо ($2$-оос бусад) нийлмэл байдаг тул бид тэгш тоог шалгахаа бүрмөсөн болиулж болно. Үүний оронд зөвхөн сондгой тоотой ажиллах хэрэгтэй.

Нэгдүгээрт, энэ нь шаардлагатай санах ойг хоёр дахин багасгах боломж олгоно. Хоёрдугаарт, алгоритмын хийх үйлдлийн тоог ойролцоогоор хоёр дахин багасгана.

Memory consumption and speed of operations¶

We should notice, that these two implementations of the Sieve of Eratosthenes use $n$ bits of memory by using the data structure vector<bool>. vector<bool> is not a regular container that stores a series of bool (as in most computer architectures a bool takes one byte of memory). It's a memory-optimization specialization of vector<T>, that only consumes $\frac{N}{8}$ bytes of memory.

Modern processors architectures work much more efficiently with bytes than with bits as they usually cannot access bits directly. So underneath the vector<bool> stores the bits in a large continuous memory, accesses the memory in blocks of a few bytes, and extracts/sets the bits with bit operations like bit masking and bit shifting.

Because of that there is a certain overhead when you read or write bits with a vector<bool>, and quite often using a vector<char> (which uses 1 byte for each entry, so 8x the amount of memory) is faster.

However, for the simple implementations of the Sieve of Eratosthenes using a vector<bool> is faster. You are limited by how fast you can load the data into the cache, and therefore using less memory gives a big advantage. A benchmark (link) shows, that using a vector<bool> is between 1.4x and 1.7x faster than using a vector<char>.

The same considerations also apply to bitset. It's also an efficient way of storing bits, similar to vector<bool>, so it takes only $\frac{N}{8}$ bytes of memory, but is a bit slower in accessing the elements. In the benchmark above bitset performs a bit worse than vector<bool>. Another drawback from bitset is that you need to know the size at compile time.

Хэсэгчилсэн шигшүүр¶

"Язгуур хүртэл шигших" оновчлолоос үзэхэд is_prime[1...n] бүхэл массивыг байнга хадгалах шаардлагагүй болно. Шигшихийн тулд зөвхөн $n$-ийн язгуур хүртэлх анхны тоонуудыг, өөрөөр хэлбэл prime[1... sqrt(n)]-г хадгалж, бүх мужийг блокуудад хувааж, блок бүрийг тусад нь шигшихэд хангалттай.

$s$ нь блокийн хэмжээг тодорхойлох тогтмол байг, тэгвэл нийт $\lceil {\frac n s} \rceil$ блок байх бөгөөд $k$ дугаар блок ($k = 0 ... \lfloor {\frac n s} \rfloor$) нь $[ks; ks + s - 1]$ хэрчим дэх тоонуудыг агуулна. Бид блокуудтай ээлжлэн ажиллаж болно, өөрөөр хэлбэл блок $k$ бүрийн хувьд бүх анхны тоог ($1$-ээс $\sqrt n$ хүртэл) давтаж, тэдгээрийг ашиглан шигшинэ. Эхний тоонуудыг боловсруулахдаа стратегиа бага зэрэг өөрчлөх шаардлагатайг тэмдэглэх нь зүйтэй: нэгдүгээрт, $[1; \sqrt n]$-ийн бүх анхны тоо өөрсдийгөө устгах ёсгүй; хоёрдугаарт, $0$ ба $1$ тоог анхны бус гэж тэмдэглэх ёстой. Сүүлийн блоктой ажиллахдаа шаардлагатай сүүлийн тоо $n$ нь блокийн төгсгөлд заавал байрлахгүй гэдгийг мартаж болохгүй.

As discussed previously, the typical implementation of the Sieve of Eratosthenes is limited by the speed how fast you can load data into the CPU caches. By splitting the range of potential prime numbers $[1; n]$ into smaller blocks, we never have to keep multiple blocks in memory at the same time, and all operations are much more cache-friendlier. As we are now no longer limited by the cache speeds, we can replace the vector<bool> with a vector<char>, and gain some additional performance as the processors can handle read and writes with bytes directly and don't need to rely on bit operations for extracting individual bits. The benchmark (link) shows, that using a vector<char> is about 3x faster in this situation than using a vector<bool>. A word of caution: those numbers might differ depending on architecture, compiler, and optimization levels.

Here we have an implementation that counts the number of primes smaller than or equal to $n$ using block sieving.

int count_primes(int n) {
    const int S = 10000;

    vector<int> primes;
    int nsqrt = sqrt(n);
    vector<char> is_prime(nsqrt + 2, true);
    for (int i = 2; i <= nsqrt; i++) {
        if (is_prime[i]) {
            primes.push_back(i);
            for (int j = i * i; j <= nsqrt; j += i)
                is_prime[j] = false;
        }
    }

    int result = 0;
    vector<char> block(S);
    for (int k = 0; k * S <= n; k++) {
        fill(block.begin(), block.end(), true);
        int start = k * S;
        for (int p : primes) {
            int start_idx = (start + p - 1) / p;
            int j = max(start_idx, p) * p - start;
            for (; j < S; j += p)
                block[j] = false;
        }
        if (k == 0)
            block[0] = block[1] = false;
        for (int i = 0; i < S && start + i <= n; i++) {
            if (block[i])
                result++;
        }
    }
    return result;
}

Блокоор шигших ажиллах хугацаа нь ердийн Эратосфены шигшүүртэй ижил (блокийн хэмжээ маш бага биш л бол), гэвч шаардлагатай санах ой $O(\sqrt{n} + S)$ болж багасах бөгөөд кэшийн үр дүн илүү сайн болно. Нөгөө талаас блок ба $[1; \sqrt{n}]$-ийн анхны тооны хос бүрт хуваалт хийгдэх бөгөөд энэ нь блокийн хэмжээ бага байх тусам муу байна. Тиймээс $S$ тогтмолыг сонгохдоо тэнцвэрийг барих шаардлагатай. Бид $10^4$-ээс $10^5$ хоорондох блокийн хэмжээнд хамгийн сайн үр дүнд хүрсэн.

Мужид анхны тоо олох¶

Заримдаа бид $R$ нь маш том байж болох (жишээ нь $1e12$) бага хэмжээтэй $[L,R]$ муж (жишээ нь $R - L + 1 \approx 1e7$) дахь бүх анхны тоог олох шаардлагатай болдог.

Ийм бодлого бодохын тулд хэсэгчилсэн шигшүүрийн санааг ашиглаж болно. Бид $\sqrt R$ хүртэлх бүх анхны тоог урьдчилан үүсгээд, тэдгээр анхны тоог ашиглан $[L, R]$ хэрчим дэх бүх нийлмэл тоог тэмдэглэнэ.

vector<char> segmentedSieve(long long L, long long R) {
    // generate all primes up to sqrt(R)
    long long lim = sqrt(R);
    vector<char> mark(lim + 1, false);
    vector<long long> primes;
    for (long long i = 2; i <= lim; ++i) {
        if (!mark[i]) {
            primes.emplace_back(i);
            for (long long j = i * i; j <= lim; j += i)
                mark[j] = true;
        }
    }

    vector<char> isPrime(R - L + 1, true);
    for (long long i : primes)
        for (long long j = max(i * i, (L + i - 1) / i * i); j <= R; j += i)
            isPrime[j - L] = false;
    if (L == 1)
        isPrime[0] = false;
    return isPrime;
}

Энэ аргын time complexity нь $O((R - L + 1) \log \log (R) + \sqrt R \log \log \sqrt R)$ юм.

Мөн бид бүх анхны тоог урьдчилан үүсгэхгүй байх боломжтой:

vector<char> segmentedSieveNoPreGen(long long L, long long R) {
    vector<char> isPrime(R - L + 1, true);
    long long lim = sqrt(R);
    for (long long i = 2; i <= lim; ++i)
        for (long long j = max(i * i, (L + i - 1) / i * i); j <= R; j += i)
            isPrime[j - L] = false;
    if (L == 1)
        isPrime[0] = false;
    return isPrime;
}

Мэдээжийн хэрэг complexity нь илүү муу буюу $O((R - L + 1) \log (R) + \sqrt R)$ болно. Гэсэн хэдий ч практикт маш хурдан ажилладаг.

Шугаман хугацааны өөрчлөлт¶

Бид алгоритмыг зөвхөн шугаман time complexity-тэй болгож өөрчилж болно. Энэ аргыг Шугаман шигшүүр өгүүлэлд тайлбарласан. Гэвч энэ алгоритм өөрийн гэсэн сул талтай.