世界总决赛选手带你玩转数论 2——质因数分解和欧几里得算法

2021-06-16 16:24:42 浏览数 (1)

笔者

笔者曾获得 ICPC 2020 世界总决赛资格，ICPC 2020 亚洲区域总决赛第五名。

本次内容

本次内容主要围绕质因数分解和欧几里得算法两部分展开，主要内容如下：

质因数分解

  • 暴力
  • Pollard Rho

欧几里得算法

  • GCD
  • 拓展 GCD
  • 类 GCD

质因数分解

对于一个质数

p

在

n!

中出现的次数为

left lfloor frac{n}{p}right rfloor left lfloor frac{n}{p^2}right rfloor left lfloor frac{n}{p^3}right rfloor cdots

暴力

暴力分解一个数

n

的时间复杂度是

O(sqrt{n})

。

代码语言：javascript复制
#define LL long long
int fac[100000];
void get_fac(LL n)
{
    int num=0;
    for (LL i=2;i*i<=n;i  )
    {
        while(n%i==0) fac[num  ]=i,n/=i;
    }
    if (n>1) fac[num  ]=n;
}

Pollard Rho

分解一个合数的高效算法。

可以将问题转化为求一个合数

n

的非平凡因子（非平凡因子是指除了

1

和本身以外的最小约数）。

从生日悖论可以得到，我们在

[1,n]

内随机

O(sqrt{n})

次整数，就有大概率出现两个相同的数。

于是我们可以考虑随机一个

x_0

出来，然后有一个多项式

f(x)

，我们令

x_i=f(x_{i-1})

。

也就是通过大概

sqrt{n}

次，我们可以得到

x_iequiv x_j (mod ; n) (i < j)

，因为

f(x)

固定，所以之后就会有

x_{i 1}equiv x_{j 1} ,x_{i 2}equiv x_{j 2},cdots

，也就是说

x_i

数列是

rho

形的，也就是会入一个环，之后一直在环内循环往复。

因为得到这样的

x_i,x_j

，我们就可以得到

|x_i-x_j|equiv 0 (mod ; n)

。

考虑如何找到这样的循环节，显然把所有的

x_i

都保存下来不太现实，因为

n

可能会比较大。

有一个方式是用两个指针，一个指针每次走一步，另一个指针每次走两步，这样当后一个指针入环以后，显然最多再进行环长步数就一定会重合（从奇偶性考虑），也就是找到循环节重合的位置时间复杂度为

O(sqrt{n})

。

而我们要找的是非平凡因子，显然合数非平凡因子一定

le sqrt{n}

，我们令

p=sqrt{n}

，那么我们只需要在

mod ; p

意义下找这样的

x_i,x_j

，这样找到一个非平凡因子的期望时间复杂度是

O(n^{frac{1}{4}})

。

代码语言：javascript复制
mt19937 mt(time(0));
LL f(LL v,LL n,LL c){LL t=mul(v,v,n) c;return t<n?t:t%n;};
LL pollard_rho(LL n,LL c)
{
    LL x = uniform_int_distribution<LL>(1,n-1)(mt),y=x;
    while (1){
        x=f(x,n,c);y=f(f(y,n,c),n,c);
        if (x==y) return n;
        LL d=gcd(abs(x-y),n);
        if (d!=1) return d;
    }
}
LL fac[100], fcnt;
void get_fac(LL n,LL cc=19260817)
{
    if (n==4){fac[fcnt  ]=2;fac[fcnt  ]=2;return;}
    if (primeQ(n)){fac[fcnt  ]=n;return;}
    LL p=n;
    while (p==n) p=pollard_rho(n,--cc);
    get_fac(p);get_fac(n/p);
}
void go_fac(LL n){fcnt=0;if (n>1) get_fac(n);}

欧几里得算法

GCD

r|a,r|b Leftrightarrow r|(a-b),r|b

于是有

(a,b)=(a-b,b)=(a; mod; b,b)
代码语言：javascript复制
void GCD(LL x,LL y)
{
    LL z=x%y;
    while(z) x=y,y=z,z=x%y;
    return y;
}

LCM

[a,b](a,b)=acdot b

需要注意的是，为了防止溢出，需要先做除法再做乘法

代码语言：javascript复制
void LCM(LL x,LL y)
{
    return x/GCD(x,y)*y;
}

裴蜀定理

对于整数

a,b,c

，如果

ax by=c

有整数解，那么一定有

(a,b)|c

拓展 GCD

用于求解

ax by=(a,b)

的任意整数解。

考虑在求 GCD 的过程中同时求出

x,y

，我们令

a'=b,b'=a; mod; b

（操作一步以后）

如果递归下去，则相当于现在已经求得

a'x b'y=(a',b')

的

x,y

了

b'=a; mod ;b =a - a/b*b

，则

bx (a-a/b*b)y=(a',b')=(a,b)Rightarrow ay b(x-a/b*y)=(a,b)

于是可以得到新的

x=y,y=x-a/b*y
代码语言：javascript复制
LL ex_gcd(LL a,LL b,LL &x,LL &y)
{
    if (b==0) {x=1;y=0;return a;}
    LL ret=ex_gcd(b,a%b,y,x);
    y-=a/b*x;
    return ret;
}

这样可以求得一个任意解

x_0,y_0

，而一般题目会对解的范围有一定的要求，可以发现，经过下面调整，仍然是一个题目的解：

x_o frac{b}{(a,b)},y_0-frac{a}{(a,b)}

。

类 GCD

类欧几里得最基本的就是求

sum_{x=0}^n left lfloor frac{ax b}{c} right rfloor

，其中

0 < a,b,c,n

，为了方便，我们用

f(a,b,c,n)

来表示这个式子。

可以先理解一下这个式子的含义，如果把

frac{ax b}{c}

看作坐标系内的一条直线方程，那么这个式子求的就是从

0le xle n

内在

x

轴上方，直线下方的整点个数。

考虑

age c

的时候，我们可以将式子分成两部分（整除的一部分可以提取出来）：

f(a,b,c,n)=frac{n(n 1)}{2}times left lfloor frac{a}{c} right rfloor (n 1) times left lfloor frac{b}{c} right rfloor f(a; mod ; c,b; mod ; c,c,n)

.

这样一来

a,b < c

了，之后根据其几何意义，我们可以考虑将直线反转到

y

轴上，变成求

y

轴到，如图所示（草图，仅为方便理解）：

这样经过反转以后，本来直线的斜率为

frac{a}{c}

，此时变成了

frac{c}{a}

，现在要求的整点个数，我们可以考虑从整个矩形减去

S

部分的整点。

为了方便表示，我们令

m= left lfloor frac{an b}{c} right rfloor

则有

因为之前

a,b < c

所以此时交换了

a,c

的位置，又回到了

age c

的情况。

因为这样的形式类似于 GCD 的求解方式，所以也被称为类欧几里得。

剩下类欧几里得还有两种形式（推导方式方式类似，就不再赘述）：

海外加速 编程算法

0 人点赞