[日常摸鱼]bzoj4802 欧拉函数-PollardRho大整数分解算法

提示：扫一扫查出行
【扫一扫了解最新限行尾号】
复制提示

啊居然要特判，卡了好久QAQ

（好像Windows下的rand和Linux下的不一样？

QwQ一些东西参考了喵铃的这篇blog：http://www.cnblogs.com/meowww/p/6400841.html （业界良心）

题目

题意：输入$n$，求$phi(n)$，$n \leq 10^{18}$

随便抽的题，刚好学习一下相关的算法。

很明显朴素的根号算法时间复杂度补滋兹，线性筛更不用想了，不过这题只需要单个欧拉函数值，还是直接考虑$phi(n)=n*\prod_{i=1}^k(1-\frac{1}{p_i})$这种求法吧，根号算法的瓶颈在于没办法快速确定有哪些约数，如果我们能快速的对$n$进行质因数分解就好了。

前置技能(pre-skill)

费马小定理

若对素数模数$p$有$a^p \equiv a \pmod p$，证明挺多的所以我就略了（上面那篇blog也有）

Miller-Rabin素性测试

对于朴素的判定素数的方法复杂度也是根号的，尽管这种复杂度其实已经挺快的了，但是有时候还需要更高效率的素性测试，注意到上面的费马小定理，如果在$[1,p-1]$内随机选十几个$a$来代到上面的式子里去，也许可以？
我实际测试了一下确实基本上能够很准确的判断素性（好吧其实我写了个对拍居然一直没拍出问题来 [真是叫人质壁分离.jpg]），在随机数据下正确率还是挺高的，之前也见过有人在比赛的时候这么判素数…
但是！有个叫卡米歇尔的人发现了这样一些数，他们是合数但是不论$a$怎么取都能通过上面这种测试，最小的一个是561=3*11*17，像这样的数还有许多，他们被称作卡米歇尔数，不过这已经不是我们要讨论的重点。
Miller-Rabin素性测试：对于某个大于1的奇数$n$，令$2^{\alpha}t=n-1$，$t$为奇数，对于$n \not | a$的某个$a$如果有$a^t \not \equiv 1 \pmod n$且$\forall 0 \leq i \leq \alpha -1,a^{2^{i}t} \not \equiv n-1 \pmod p$，则$n$一定是合数。可以证明和它等价的逆否命题，具体的上面那篇blog讲的非常清楚就不债述了QwQ
这就是miller-rabin素性测试了，看起来比较麻烦而且，不过这种方法的好处就在于它不存在像上面卡米歇尔数一样不论如何都无法正确判断的情况。需要注意的是实现的时候由于经常是对很大的数进行素性测试，乘法需要用类似快速幂那样的写法来实现，不然可能会溢出，具体的实现放在后面的代码里面了，就是那个check函数。

Pollard-Rho大整数分解算法(zheng ti)

好了现在终于可以回到最开始的那个问题了，我们需要快速地对一个数进行质因数分解。

我们假装找到了一个因子

分解之前先用上面的素性测试判断$n$是不是一个素数，如果不是的话，假装我们已经可以找到了$n$的一个因子$d$，那么就可以递归的去分解$d$和$n/d$，现在我们的问题就是要怎么快速的找到这个因子。
嗯还是上面的那篇blog，已经非常清楚的讲了怎么分解出这个$d$了~我不觉得我可以写得更好所以我就不写了（才不是因为懒）

好了现在回到原来的题目，用Pollard-Rho算法分解出$n$的所有因子排个序就好辣。

注意：特判1 特判1 特判1

细节都在里面啦，我写的好像比较挫在bzoj上跑了一百多ms…

T_T搞不懂那些10ms内的怎么写的

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<queue>
#include<algorithm>
using namespace std;
typedef long long lint;inline lint mul(lint a,lint b,lint p)
{
lint res=0;a%=p;b%=p;
for(;b;b>>=1,a=(a<<1)%p)if(b&1)res=(res+a)%p;
return res;
}
inline lint pow_mod(lint a,lint b,lint p)
{
lint res=1;a%=p;
for(;b;b>>=1,a=mul(a,a,p))if(b&1)res=mul(res,a,p);
return res;
}
inline lint gcd(lint a,lint b)
{
return !b?a:gcd(b,a%b);
}
inline lint get(lint x)
{
return rand()%x;
//return mul(mul(rand(),rand(),x),mul(rand(),rand(),x),x);
}
inline lint range(lint a,lint b)
{
return a+get(b-a+1);
}
inline bool check(lint x)
{
lint temp[15]={2,3,5,7,11,13,17,19,23,29,31,37,41,43,47};
for(register int i=0;i<15;i++)
{
if(x==temp[i])return 1;
if(x%temp[i]==0)return 0;
}
lint s,a,k=0,t=x-1;
while(!(t&1))
{
t>>=1;
k++;
}
for(register lint i=1;i<=20;i++)
{
a=range(2,x-1);
s=pow_mod(a,t,x);
if(s==1)return 1;
for(register lint j=0;j<=k-1;j++)
{
s=mul(s,s,x);
if(s==x-1)return 1;
}
}
return 0;
}lint q[115],tot;
#define f(p) ((mul(p,p,x)+c)%x)
inline void solve(lint x)
{
if(check(x))
{
q[++tot]=x;return;
}
while(1)
{
lint a,b,c;
a=b=range(0,x-1);c=range(0,x-1);b=f(b);
while(a!=b)
{
lint temp=a-b;
temp=gcd(abs(temp),x);
if(1<temp&&temp<x)
{
solve(temp);solve(x/temp);
return;
}
a=f(a);b=f(f(b));
}
}
}
int main()
{
srand(19260817);
lint n,res;scanf("%lld",&n);
if(n==1)
{
printf("1");
return 0;
}
solve(n);sort(q+1,q+tot+1);
res=n;
for(register int i=1;i<=tot;i++)if(q[i]!=q[i-1])res=(res/q[i])*(q[i]-1);
printf("%lld",res);
return 0;
}