Агуулгыг алгасах

Бүхэл тоог үржигдэхүүнд задлах

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

Хэрэв та задлахыг хүсэж буй тоо чинь үнэндээ анхны тоо бол ихэнх алгоритм маш удаан ажиллана гэдгийг анхаараарай. Энэ нь ялангуяа Фермагийн, Поллардын p-1 болон Поллардын ро задаргааны алгоритмуудад үнэн. Тиймээс тоог задлахыг оролдохын өмнө магадлалын (эсвэл хурдан детерминистик) анхны тооны шалгуур хийх нь хамгийн утга учиртай.

Туршилтын хуваалт

Энэ бол анхны үржигдэхүүнд задаргааг олох хамгийн энгийн алгоритм юм.

Бид боломжит хуваагч $d$ бүрээр хуваана. Нийлмэл тоо $n$-ийн бүх анхны үржигдэхүүн $\sqrt{n}$-ээс их байх боломжгүйг ажиглаж болно. Тиймээс бид зөвхөн $2 \le d \le \sqrt{n}$ хуваагчдыг шалгах хэрэгтэй бөгөөд энэ нь бидэнд анхны үржигдэхүүнд задаргааг $O(\sqrt{n})$-д өгнө. (Энэ бол псевдо-полиномын хугацаа, өөрөөр хэлбэл оролтын утгын хувьд полиномын боловч оролтын битийн тооны хувьд экспоненциал.)

Хамгийн бага хуваагч нь анхны тоо байх ёстой. Бид задарсан тоог хасаад үйл ажиллагааг үргэлжлүүлнэ. Хэрэв $[2; \sqrt{n}]$ мужид ямар ч хуваагч олохгүй бол тоо өөрөө анхны байх ёстой.

vector<long long> trial_division1(long long n) {
    vector<long long> factorization;
    for (long long d = 2; d * d <= n; d++) {
        while (n % d == 0) {
            factorization.push_back(d);
            n /= d;
        }
    }
    if (n > 1)
        factorization.push_back(n);
    return factorization;
}

Дугуйт задаргаа

Энэ бол туршилтын хуваалтын оновчлол юм. Тоо 2-т хуваагдахгүй гэдгийг мэдсэн бол бусад тэгш тоог шалгах шаардлагагүй. Ингэснээр бидэнд шалгах тооны зөвхөн $50\%$ үлдэнэ. 2-ыг ялган авч сондгой тоо гаргасны дараа бид зүгээр л 3-аас эхэлж, зөвхөн бусад сондгой тоог тоолж болно.

vector<long long> trial_division2(long long n) {
    vector<long long> factorization;
    while (n % 2 == 0) {
        factorization.push_back(2);
        n /= 2;
    }
    for (long long d = 3; d * d <= n; d += 2) {
        while (n % d == 0) {
            factorization.push_back(d);
            n /= d;
        }
    }
    if (n > 1)
        factorization.push_back(n);
    return factorization;
}

Энэ аргыг цааш өргөтгөж болно. Хэрэв тоо 3-т хуваагдахгүй бол дараагийн тооцоололд 3-ын бусад бүх үржвэрийг мөн үл тоомсорлож болно. Тиймээс бид зөвхөн $5, 7, 11, 13, 17, 19, 23, \dots$ тоонуудыг шалгах хэрэгтэй. Эдгээр үлдсэн тоонуудын хэв маягийг ажиглаж болно. Бид $d \bmod 6 = 1$ ба $d \bmod 6 = 5$ байх бүх тоог шалгах хэрэгтэй. Ингэснээр бидэнд шалгах тооны ердөө $33.3\%$ хувь үлдэнэ. Үүнийг эхлээд 2 ба 3 анхны тоог ялган авч, дараа нь 5-аас эхлэн $6$ модулиар $1$ ба $5$ үлдэгдлийг л тоолох замаар хэрэгжүүлж болно.

Here is an implementation for the prime number 2, 3 and 5. It is convenient to store the skipping strides in an array.

vector<long long> trial_division3(long long n) {
    vector<long long> factorization;
    for (int d : {2, 3, 5}) {
        while (n % d == 0) {
            factorization.push_back(d);
            n /= d;
        }
    }
    static array<int, 8> increments = {4, 2, 4, 2, 4, 6, 2, 6};
    int i = 0;
    for (long long d = 7; d * d <= n; d += increments[i++]) {
        while (n % d == 0) {
            factorization.push_back(d);
            n /= d;
        }
        if (i == 8)
            i = 0;
    }
    if (n > 1)
        factorization.push_back(n);
    return factorization;
}

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

Урьдчилан тооцоолсон анхны тоо

Дугуйт задаргааны аргыг хязгааргүй өргөтгөвөл бидэнд шалгах зөвхөн анхны тоонууд үлдэнэ. Үүнийг шалгах сайн арга бол Эратосфены шигшүүр-ээр $\sqrt{n}$ хүртэлх бүх анхны тоог урьдчилан тооцоолж, тэдгээрийг тус тусад нь шалгах явдал юм.

vector<long long> primes;

vector<long long> trial_division4(long long n) {
    vector<long long> factorization;
    for (long long d : primes) {
        if (d * d > n)
            break;
        while (n % d == 0) {
            factorization.push_back(d);
            n /= d;
        }
    }
    if (n > 1)
        factorization.push_back(n);
    return factorization;
}

Фермагийн задаргааны арга

Сондгой нийлмэл тоо $n = p \cdot q$-г хоёр квадратын ялгавар $n = a^2 - b^2$ хэлбэрээр бичиж болно:

$$n = \left(\frac{p + q}{2}\right)^2 - \left(\frac{p - q}{2}\right)^2$$

Фермагийн задаргааны арга нь эхний квадрат $a^2$-г таамаглаж, үлдсэн хэсэг $b^2 = a^2 - n$ мөн квадрат тоо эсэхийг шалгах замаар энэ баримтыг ашиглахыг оролддог. Хэрэв тийм бол бид $n$-ийн $a - b$ ба $a + b$ үржигдэхүүнүүдийг олсон болно.

int fermat(int n) {
    int a = ceil(sqrt(n));
    int b2 = a*a - n;
    int b = round(sqrt(b2));
    while (b * b != b2) {
        a = a + 1;
        b2 = a*a - n;
        b = round(sqrt(b2));
    }
    return a - b;
}

Хэрэв $p$ ба $q$ хоёр үржигдэхүүний ялгаа бага бол энэ задаргааны арга маш хурдан байж болно. Алгоритм $O(|p - q|)$ хугацаанд ажиллана. Гэвч практикт энэ аргыг ховор ашигладаг. Үржигдэхүүнүүд бие биенээсээ хол болмогц энэ нь маш удаан болно.

Гэсэн хэдий ч энэ аргад хамаарах олон оновчлолын сонголт байсаар байна. Тогтмол бага тоогоор модуль авсан $a^2$ квадратуудыг харснаар зарим $a$ утгуудыг үзэх шаардлагагүйг ажиглаж болно, учир нь тэдгээр нь $a^2 - n$ квадрат тоо үүсгэж чадахгүй.

Поллардын $p - 1$ арга

$n$ тоо бага $\mathrm{B}$-ийн хувьд $p - 1$ нь $\mathrm{B}$-зэрэг гөлгөр байх дор хаяж нэг анхны үржигдэхүүн $p$-тэй байх магадлал маш өндөр. Бүхэл тоо $m$$\mathrm{B}$-зэрэг гөлгөр гэж хэлнэ, хэрэв $m$-г хуваадаг анхны зэрэг бүр хамгийн ихдээ $\mathrm{B}$ байвал. Албан ёсоор $\mathrm{B} \geqslant 1$ байг, мөн $m$ нь дурын эерэг бүхэл тоо байг. $m$-ийн анхны үржигдэхүүнд задаргаа нь $m = \prod {q_i}^{e_i}$ байг, энд $q_i$ бүр анхны тоо, $e_i \geqslant 1$. Тэгвэл бүх $i$-гийн хувьд ${q_i}^{e_i} \leqslant \mathrm{B}$ байвал $m$ нь $\mathrm{B}$-зэрэг гөлгөр болно. Жишээ нь $4817191$-ийн анхны үржигдэхүүнд задаргаа нь $1303 \cdot 3697$. Мөн $1303 - 1$ ба $3697 - 1$ утгууд нь тус тус $31$-зэрэг гөлгөр ба $16$-зэрэг гөлгөр юм, учир нь $1303 - 1 = 2 \cdot 3 \cdot 7 \cdot 31$ ба $3697 - 1 = 2^4 \cdot 3 \cdot 7 \cdot 11$. 1974 онд Жон Поллард нийлмэл тооноос $p-1$ нь $\mathrm{B}$-зэрэг гөлгөр байх $p$ үржигдэхүүнийг гаргаж авах аргыг зохион бүтээсэн.

Санаа нь Фермагийн бага теорем-оос гарна. $n$-ийн задаргаа нь $n = p \cdot q$ байг. Хэрэв $a$ нь $p$-тэй харилцан анхны бол дараах мэдэгдэл биелнэ гэж хэлдэг:

$$a^{p - 1} \equiv 1 \pmod{p}$$

Энэ нь мөн дараахыг илэрхийлнэ

$${\left(a^{(p - 1)}\right)}^k \equiv a^{k \cdot (p - 1)} \equiv 1 \pmod{p}.$$

Тиймээс $p - 1 ~|~ M$ байх дурын $M$-ийн хувьд $a^M \equiv 1$ болохыг бид мэднэ. Энэ нь $a^M - 1 = p \cdot r$ гэсэн үг бөгөөд иймээс $p ~|~ \gcd(a^M - 1, n)$ болно.

Тиймээс $n$-ийн үржигдэхүүн $p$-ийн хувьд $p - 1$ нь $M$-г хуваадаг бол бид Евклидийн алгоритм ашиглан үржигдэхүүнийг гаргаж авч болно.

$\mathrm{B}$-зэрэг гөлгөр тоо бүрийн үржвэр болох хамгийн бага $M$ нь $\text{lcm}(1,~2~,3~,4~,~\dots,~B)$ болох нь тодорхой. Эсвэл өөрөөр:

$$M = \prod_{\text{prime } q \le B} q^{\lfloor \log_q B \rfloor}$$

Хэрэв $n$-ийн бүх анхны үржигдэхүүн $p$-ийн хувьд $p-1$ нь $M$-г хуваадаг бол $\gcd(a^M - 1, n)$ зүгээр л $n$ болохыг анхаараарай. Энэ тохиолдолд бид үржигдэхүүн авахгүй. Тиймээс бид $M$-г тооцоолж байх зуураа $\gcd$-г олон удаа хийхийг оролдоно.

Зарим нийлмэл тоо бага $\mathrm{B}$-ийн хувьд $p-1$ нь $\mathrm{B}$-зэрэг гөлгөр байх $p$ үржигдэхүүнгүй байдаг. Жишээ нь $100~000~000~000~000~493 = 763~013 \cdot 131~059~365~961$ нийлмэл тооны хувьд $p-1$ утгууд нь харгалзан $190~753$-зэрэг гөлгөр ба $1~092~161~383$-зэрэг гөлгөр байна. Энэ тоог задлахын тулд бид $B \geq 190~753$ сонгох хэрэгтэй болно.

In the following implementation we start with $\mathrm{B} = 10$ and increase $\mathrm{B}$ after each each iteration.

long long pollards_p_minus_1(long long n) {
    int B = 10;
    long long g = 1;
    while (B <= 1000000 && g < n) {
        long long a = 2 + rand() %  (n - 3);
        g = gcd(a, n);
        if (g > 1)
            return g;

        // compute a^M
        for (int p : primes) {
            if (p >= B)
                continue;
            long long p_power = 1;
            while (p_power * p <= B)
                p_power *= p;
            a = power(a, p_power, n);

            g = gcd(a - 1, n);
            if (g > 1 && g < n)
                return g;
        }
        B *= 2;
    }
    return 1;
}

Энэ бол магадлалын алгоритм гэдгийг анхаараарай. Үүний үр дагавар нь алгоритм ямар ч үржигдэхүүн олж чадахгүй байх боломж байдагт оршино.

Хүндрэл нь давталт тутамд $O(B \log B \log^2 n)$ юм.

Поллардын ро алгоритм

Поллардын Ро алгоритм бол Жон Поллардын бас нэгэн задаргааны алгоритм юм.

Тооны анхны үржигдэхүүнд задаргаа нь $n = p q$ байг. Алгоритм нь $\{x_i\} = \{x_0,~f(x_0),~f(f(x_0)),~\dots\}$ псевдо-санамсаргүй дарааллыг авч үздэг, энд $f$ нь олон гишүүнт функц, ихэвчлэн $c = 1$-тэйгээр $f(x) = (x^2 + c) \bmod n$ сонгодог.

Энэ тохиолдолд бид $\{x_i\}$ дараалалд сонирхолгүй. Бид $\{x_i \bmod p\}$ дараалалд илүү сонирхолтой. $f$ нь олон гишүүнт функц бөгөөд бүх утга $[0;~p)$ мужид байдаг тул энэ дараалал эцэст нь цикл рүү нийлнэ. Төрсөн өдрийн парадокс нь давталт эхлэх хүртэл хүлээгдэж буй элементийн тоо $O(\sqrt{p})$ байхыг харуулдаг. Хэрэв $p$ нь $\sqrt{n}$-ээс бага бол давталт $O(\sqrt[4]{n})$-д эхлэх магадлалтай.

$n = 2206637$, $p = 317$, $x_0 = 2$ ба $f(x) = x^2 + 1$-тэй ийм $\{x_i \bmod p\}$ дарааллын дүрслэлийг энд үзүүлэв. Дарааллын хэлбэрээс алгоритмыг яагаад Поллардын $\rho$ алгоритм гэж нэрлэдгийг маш тодорхой харж болно.

Поллардын ро алгоритмын дүрслэл

Гэвч нэг нээлттэй асуулт байсаар байна. $p$ тоог өөрийг нь мэдэхгүйгээр $\{x_i \bmod p\}$ дарааллын шинж чанарыг хэрхэн ашиглах вэ?

Үнэндээ энэ нь нэлээд хялбар. $\{x_i \bmod p\}_{i \le j}$ дараалалд цикл байх зайлшгүй бөгөөд хүрэлцээтэй нөхцөл нь $x_s \equiv x_t \bmod p$ байх $s, t \le j$ хоёр индекс байх явдал юм. Энэ тэгшитгэлийг $x_s - x_t \equiv 0 \bmod p$ гэж дахин бичиж болох бөгөөд энэ нь $p ~|~ \gcd(x_s - x_t, n)$-тэй адил юм.

Тиймээс хэрэв бид $g = \gcd(x_s - x_t, n) > 1$ байх $s$ ба $t$ хоёр индекс олвол цикл болон $n$-ийн $g$ үржигдэхүүнийг олсон болно. $g = n$ байх боломжтой. Энэ тохиолдолд бид жинхэнэ үржигдэхүүн олоогүй тул алгоритмыг өөр параметртэйгээр (өөр эхлэх утга $x_0$, олон гишүүнт функц $f$ дэх өөр тогтмол $c$) давтах ёстой.

Циклийг олохын тулд бид ямар ч түгээмэл цикл илрүүлэх алгоритм ашиглаж болно.

Флойдын цикл олох алгоритм

Энэ алгоритм нь дараалал дээгүүр өөр өөр хурдтай хөдөлдөг хоёр заагч ашиглан циклийг олдог. Давталт бүрийн үед эхний заагч нэг элемент урагшлах бол хоёр дахь заагч нэг алгасаж элемент бүр рүү урагшилна. Энэ санааг ашиглан хэрэв цикл байвал ямар нэг үед хоёр дахь заагч давталтын явцад эхнийхтэй уулзана гэдгийг харахад амархан. Хэрэв циклийн урт нь $\lambda$ бөгөөд $\mu$ нь цикл эхлэх эхний индекс бол алгоритм $O(\lambda + \mu)$ хугацаанд ажиллана.

Энэ алгоритмыг мөн яст мэлхий (удаан заагч) ба туулай (хурдан заагч) хоёр уралдсан үлгэрт үндэслэн Яст мэлхий ба туулайн алгоритм гэж нэрлэдэг.

Үнэндээ энэ алгоритмыг ашиглан $\lambda$ ба $\mu$ параметрүүдийг тодорхойлох боломжтой (мөн $O(\lambda + \mu)$ хугацаа ба $O(1)$ санах ойд). Цикл илэрсэн үед алгоритм 'True' буцаана. Хэрэв дараалалд цикл байхгүй бол функц төгсгөлгүй давтагдана. Гэвч Поллардын Ро алгоритмыг ашиглан үүнээс сэргийлж болно.

function floyd(f, x0):
    tortoise = x0
    hare = f(x0)
    while tortoise != hare:
        tortoise = f(tortoise)
        hare = f(f(hare))
    return true

Implementation

First, here is an implementation using the Floyd's cycle-finding algorithm. The algorithm generally runs in $O(\sqrt[4]{n} \log(n))$ time.

long long mult(long long a, long long b, long long mod) {
    return (__int128)a * b % mod;
}

long long f(long long x, long long c, long long mod) {
    return (mult(x, x, mod) + c) % mod;
}

long long rho(long long n, long long x0=2, long long c=1) {
    long long x = x0;
    long long y = x0;
    long long g = 1;
    while (g == 1) {
        x = f(x, c, n);
        y = f(y, c, n);
        y = f(y, c, n);
        g = gcd(abs(x - y), n);
    }
    return g;
}

The following table shows the values of $x$ and $y$ during the algorithm for $n = 2206637$, $x_0 = 2$ and $c = 1$.

$$ \newcommand\T{\Rule{0pt}{1em}{.3em}} \begin{array}{|l|l|l|l|l|l|} \hline i & x_i \bmod n & x_{2i} \bmod n & x_i \bmod 317 & x_{2i} \bmod 317 & \gcd(x_i - x_{2i}, n) \\ \hline 0 & 2 & 2 & 2 & 2 & - \\ 1 & 5 & 26 & 5 & 26 & 1 \\ 2 & 26 & 458330 & 26 & 265 & 1 \\ 3 & 677 & 1671573 & 43 & 32 & 1 \\ 4 & 458330 & 641379 & 265 & 88 & 1 \\ 5 & 1166412 & 351937 & 169 & 67 & 1 \\ 6 & 1671573 & 1264682 & 32 & 169 & 1 \\ 7 & 2193080 & 2088470 & 74 & 74 & 317 \\ \hline \end{array}$$

The implementation uses a function mult, that multiplies two integers $\le 10^{18}$ without overflow by using a GCC's type __int128 for 128-bit integer. If GCC is not available, you can using a similar idea as binary exponentiation.

long long mult(long long a, long long b, long long mod) {
    long long result = 0;
    while (b) {
        if (b & 1)
            result = (result + a) % mod;
        a = (a + a) % mod;
        b >>= 1;
    }
    return result;
}

Alternatively you can also implement the Montgomery multiplication.

As stated previously, if $n$ is composite and the algorithm returns $n$ as factor, you have to repeat the procedure with different parameters $x_0$ and $c$. E.g. the choice $x_0 = c = 1$ will not factor $25 = 5 \cdot 5$. The algorithm will return $25$. However, the choice $x_0 = 1$, $c = 2$ will factor it.

Брентийн алгоритм

Брент Флойдтой төстэй аргыг хоёр заагч ашиглан хэрэгжүүлдэг. Ялгаа нь заагчдыг харгалзан нэг ба хоёр байрлалаар урагшлуулахын оронд 2-ын зэргээр урагшлуулдагт оршино. $2^i$ нь $\lambda$ ба $\mu$-ээс их болмогц бид циклийг олно.

function floyd(f, x0):
    tortoise = x0
    hare = f(x0)
    l = 1
    while tortoise != hare:
        tortoise = hare
        repeat l times:
            hare = f(hare)
            if tortoise == hare:
                return true
        l *= 2
    return true

Брентийн алгоритм мөн шугаман хугацаанд ажилладаг боловч $f$ функцийн тооцооллыг бага ашигладаг тул ерөнхийдөө Флойдынхоос хурдан байдаг.

Implementation

The straightforward implementation of Brent's algorithm can be sped up by omitting the terms $x_l - x_k$ if $k < \frac{3 \cdot l}{2}$. In addition, instead of performing the $\gcd$ computation at every step, we multiply the terms and only actually check $\gcd$ every few steps and backtrack if overshot.

long long brent(long long n, long long x0=2, long long c=1) {
    long long x = x0;
    long long g = 1;
    long long q = 1;
    long long xs, y;

    int m = 128;
    int l = 1;
    while (g == 1) {
        y = x;
        for (int i = 1; i < l; i++)
            x = f(x, c, n);
        int k = 0;
        while (k < l && g == 1) {
            xs = x;
            for (int i = 0; i < m && i < l - k; i++) {
                x = f(x, c, n);
                q = mult(q, abs(y - x), n);
            }
            g = gcd(q, n);
            k += m;
        }
        l *= 2;
    }
    if (g == n) {
        do {
            xs = f(xs, c, n);
            g = gcd(abs(xs - y), n);
        } while (g == 1);
    }
    return g;
}

Бага анхны тоонуудын туршилтын хуваалтыг Поллардын ро алгоритмын Брентийн хувилбартай хослуулснаар маш хүчирхэг задаргааны алгоритм болно.

Дасгал бодлогууд