整数分解及素性判断学习笔记

整数分解

    整数分解事实上有三个子问题：

    1. 合性判断。即判断一个数是否为合数。

    2. 素性判断。即确切判断一个数是否为素数。

    3. 寻找可能的质因子。  

一个大概率的素性判断——Miller-Rabin算法

费马小定理

    若p为质数，那么。

    当我们要判断整数p是否为素数时，若，则p为合数。

    但若找不到满足条件的x，p只有99.999745%的概率为素数。因为存在一种数——Carmichael卡迈克尔数。

卡迈克尔数

    若p为合数，且，则称p为卡迈克尔数。

    卡迈克尔数有至少三个不同素因子。

    1e8以内只有255个卡迈克尔数，前三个为561，1105，1729。

    此时考虑一个新的定理：欧几里得引理。

欧几里得引理

    若p为素数，且,，则。

    证明：

       

    首先将p写成的形式，那么对于已选择的x，首先运算，然后将这个数平方，得，重复此步骤k-1次，得，即。

    那么在这k步中，若出现了1，而上一次的结果不是1或-1，则p为合数；若最后结果不为1，则p为合数。

    那么当p通过一次x的测试，则p为合数的概率为1/4；那么若p通过5次测试，素性判断正确的概率为99.9%。

    如果选择x=2,3,5,7，那么在2.5*10^13内仅有3215031751会判断错误。

一个复杂度并不太好的素性判断

    除了2，所有素数都为奇数。对于一个奇数p，若它为合数，则可以写成的形式。令，则。若为完全平方数，则p为合数。

    但这个算法十分不优秀..

优化掉这个不太优秀的算法——Pollard Rho

    考虑这样一个问题：有没有办法快速且高概率的找到一个小因子呢。

    很容易发现一个不高速低概率的方式：随机试除。

    这个算法概率极低，基本相当于从[1,N]随机选择一个数，这个数=x的概率——1/N。

    但考虑这样一个算法：从[1,N]随机选择两个数a和b，|a-b|=x的概率。显然是2/N。

    这个优化非常小，但我们再优化一下这个算法：从[1,N]随机选择k个数x1...xk，存在一对i和j使得|xi-xj|=x的概率。

    此时考虑一个经典问题——Birthday paradox生日悖论。

生日悖论

    随机选择N个人，设事件P为这N个人中存在两个人生日相同。求P发生的概率。

    考虑没有两人生日相同的概率。将人一个一个加入，第二个人和第一个人生日不同的概率为，第三个人和前两个人生日不同的概率为...，第n个人和前面所有人生日不同的概率为。

   

    当时，P发生的概率大于50%。

    考虑推广生日悖论。生日悖论等价于N个数中随机选k个数，存在一对数差为0的概率。可以将其推广为：N个数中随机选k个数，存在一对数差为x的概率。即上方的问题。

    也就是说，如果我们选择个数，就可以50%的概率找到N的一个小因子。但是可能太过巨大，无法存在内存中，需要特殊处理。

    考虑一种流水线的随机方式和floyd判环方法。

Pollard推荐的随机数生成方法

    定义。a可以随机生成或直接用a=2。那么随机选取后，，每次和前一个值作差即可。

    但当伪随机数出现循环节，后面的伪随机数就没有意义了，我们需要判断掉这种情况。

    这里用到的判断方法为floyd算法。

    考虑跑道上的两个人，一个人的速度为另一个人的两倍。那么当跑道是环形，第二个人一定可以追上第一个人。这样就可以不消耗空间的判断是否出现了循环节。

int find_factorplus(int N) {
    a = 2;
    b = a;
    do {
        a = f(a);//a runs once
        b = f(f(b));//b runs twice as fast
        p = GCD( abs( b - a ) , N);
        if( p > 1 ) return p;//Found factor: p
    } while( b != a );
    return 0;//Failed. :-(
}

    这样我们就得到了Pollard Rhd算法。

Pollard Rho与Miller-Rabin的结合——一个优秀的整数分解算法

    Pollard Rho算法的优势是可以在的复杂度内得到一个因子p；而当N的因子少而大时，我们可以使用M_R算法来验证下去。这样就得到了一个大概率且高效的整数分解算法。

#include <bits/stdc++.h>
using namespace std;

typedef long long LL;

const int ps[] = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37};
const int pcnt = sizeof(ps) / sizeof(int);

inline LL quick_mul(LL a, LL b, LL m)
{
    LL d = ((long double)a / m * b + 1e-8);
    LL r = a * b - d * m;
    return r < 0 ? r + m : r;
}

inline LL quick_pow(LL a,LL b,LL m)
{
    LL res = 1;
    for(; b; b >>= 1, a = quick_mul(a, a, m))
        if (b & 1)
            res = quick_mul(res, a, m);
    return res;
}

inline LL gcd(LL a, LL b)
{
    if (!a || !b)
        return a + b;
    int t = __builtin_ctzll(a | b);
    a >>= __builtin_ctzll(a);
    do
    {
        b >>= __builtin_ctzll(b);
        if (a > b)
        {
            LL t = b;
            b = a;
            a = t;
        }
        b -= a;
    }
    while (b != 0);
    return a << t;
}

inline int Miller_Rabin(LL n)
{
    if (n == 1)
        return 0;
    if (n == 2 || n == 3 || n == 5)
        return 1;
    if (!(n & 1) || (n % 3 == 0) || (n % 5 == 0))
        return 0;
    LL m = n - 1;
    int k = 0;
    while (!(m & 1))
        m >>= 1, ++k;
    for (int i = 0; i < pcnt && ps[i] < n; ++i)
    {
        LL x = quick_pow(ps[i], m, n), y = x;
        for (int i = 0; i < k; ++i)
        {
            x = quick_mul(x, x, n);
            if (x == 1 && y != 1 && y != n - 1)
                return 0;
            y = x;
        }
        if (x != 1)
            return 0;
    }
    return 1;
}

inline LL next_rand(LL x, LL n, LL a)
{
    LL t = quick_mul(x, x, n) + a;
    return t < n ? t : t - n;
}

const int M = (1 << 7) - 1;
inline LL Pollard_Rho(LL n)
{
    if (n % 2 == 0)
        return 2;
    if (n % 3 == 0)
        return 3;
    LL x = 0, y = 0, t = 1, q = 1, a = (rand() % (n - 1)) + 1;
    for (int k = 2; true; k <<= 1, y = x, q = 1)
    {
        for (int i = 1; i < k; ++i)
        {
            x = next_rand(x, n, a);
            q = quick_mul(q, abs(x - y), n);
            if (!(i & M))
            {
                t = gcd(q, n);
                if (t > 1)
                    break;
            }
        }
        if (t > 1 || (t = gcd(q, n)) > 1)
            break;
    }
    if (t == n)
        for (t = 1; t == 1; t = gcd(abs((x = next_rand(x, n, a)) - y), n));
    return t;
}

LL f[105];
int cnt;

void solve(LL n)
{
    if (n == 1)
        return;
    if (Miller_Rabin(n))
    {
        f[cnt++] = n;
        return;
    }
    LL t = n;
    while (t == n)
        t = Pollard_Rho(n);
    solve(t);
    solve(n / t);
}

int main()
{
    LL n;
    while (~scanf("%lld", &n))
    {
        cnt = 0;
        solve(n);
        sort(f, f + cnt);
        for (int i = 0; i < cnt; ++i)
            printf("%lld ", f[i]);
        putchar('
');
    }
    return 0;
}

原根的存在性及个数证明

   https://www.cnblogs.com/zhixingr/p/6822370.html