【JZOJ5745】幂

Time Limits: 3000 ms Memory Limits: 1048576 KB Detailed Limits
Description

Input

Output

Sample

Input1:
3
Output1:
1
Input2:
4
Output2:
-1
Input3:
15
Output3:
2

Data Constraint

---

analysis

$$\forall a^{n^k}\equiv a(\mod n)$$
很容易发现当n有非1的平方因子时一定无解，因为设 $n=p^2q(p>1)$，当 $a=p$时，又因为 $n\geq2$，那么肯定找不到k
那么有 $n=\prod p_i$
该式子可以分解为多个 $a^{n^k}\equiv a(\mod p_i)$
则会有 $n^k\equiv 1(\mod p_i-1)$
因此，设 $L=LCM(p_1-1,p_2-1,\cdots ,p_s-1)$
那么 $$n^k\equiv 1(\mod L)$$
同时当 $\gcd(n,L)>1$也无解

又$\because n^{\varphi(L)}\equiv 1(\mod L)$

$$k|\varphi(L)$$
枚举 $\varphi(L)$的因数即可
可能用到快速幂和快速加。

等等。。。您1e18的值域要找约数？？
其实只需要将其质因数分解再暴力约数即可。

如何快速质因数分解呢？

这里需要用到$Miller\_Rabin$和$Pollard\_rho$

$Miller\_Rabin$快速判定质数法
说到质数，很容易想到费马小定理$a^{p-1}\equiv 1 (\mod p)\ \ \ \ \ (1\leq a\lt p)$
难道用这个判就好了吗？
很可惜不能，有一些强伪素数对于任意的a都满足这个式子，但是他是个合数，如561
但是有一个判定质数的充分必要条件，那就是$\forall_{1\lt a\lt p-1}a^2\not\equiv 1(\mod p)$
那么，令$p-1=2^uv$
只用随机a，判断$a^p-1\equiv 1(\mod p)$且
$a^{2^{u'}v}\equiv 1(\mod p)$同时$a^{2^{u'}v}\equiv \pm1(\mod p)$
一次随机的正确率为$\frac 34$则随机个10遍就好了。

bool check(ll x,ll u,int t,ll mo){
    ll n=power(x,u);//power对mo取模
    if(n==1 || n==mo-1)return 1;
    for(int i=1;i<=t;i++){
        if(n==mo-1)return 1;
        n=n*n%mo;
        if(n==1)return 0;
    }return 0;
}
bool miller_rabin(ll x){
    if(x<2)return 0;
    if(x==2)return 1;
    if(x%2==0)return 0;
    ll u=x-1,k;int s=0;
    while(u%2==0)u/=2,++s;
    for(int i=1;i<=10;i++){
        do k=rand()%x;while(k==0);
        if(!check(k,u,s,x))return 0;
    }return 1;
}

$Pollard\_rho$快速找因数
首先先引入生日悖论
若要从1~1000中选取1，求选取成功的几率？
答案是$\frac 1{1000}$
但是如果变为随机两个数，求其差为1的概率为多少，这就变成了$\frac 1{500}$
如果取三个数，存在两者的差为1的概率，将增大为$\frac 3{500}$
通过实验，若选30个数概率将超过$50\%$
重复实验证明，当1~n中选择k个数，使得存在两个数差值为某值，若$k\leq\sqrt n$则其概率将大于$50\%$

其次，如何将生日悖论引入找因数当中呢？
假设先随便弄出一个数列，判断两两差的绝对值是否为n的因数？
概率较小了一点
因为我们知道$gcd(x,n)|n$所以我们可以通过这种方式去找因数
弄出k个数，对两两的差与n取gcd，若值不为1或n则求出的gcd一定是它的平凡的因数
难道这样就可以吗？能够发现，这样做选取的k可以为$n^{\frac 14}$但空间有一点大

Pollard’s rho算法
其实并不需要将k个数弄出来，只需一个一个地生成并检查连续的两个数，反复执行这个步骤并希望能够得到我们想要的数。
能够符合这个条件的变换则是伪随机$f(x)=(x^2+c)\mod n$
也就是$x_{i+1}=f(x_i)$拿$f_x$和另外一个x进行判断，而另外一个也要求能通过之前的变换而来
这种方案最终是存在循环的，如何判断循环呢
要个很巧妙的方法
这里拿$x_{2^i+1\to2^{i+1}-1}$和$x_{2^{i+1}}$进行判断，若判断的过程中，存在$x_{2^i+k}=x_{2^{i+1}}$，那么说明已经至少走了一圈

ll pollard_rho(ll x,ll c){
    int i=1,k=2;
    ll x0=rand()%x,y=x0;
    for(;;){
        ++i;
        x0=(x0*x0+c)%x;
        ll d=gcd(abs(x0-y),x);
        if(d!=1 && d!=x)return d;
        if(y==x0)return pollard_rho(x,rand()%x+1);//x0错误
        if(i==k)k<<=1,y=x0;
    }
}

这样一来，分解质因数的过程就变为：
1. 判断n>1
2. 判断n是否为质数(Miller Rabin)
3. 找一个因数d
4. 分成两部分d和n/d继续分解

#include<cstring>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<ctime>
#include<cstdlib>
#define N 100100
#define ll long long

using namespace std;

int t,ap[N];
ll p[N],n,ans,lcm,phi;

ll R(){
    rand();
    ll s=0;for(int i=1;i<=5;i++)s=s*100+rand()%100;
    return s+13;
}
ll qmul(ll a,ll b,ll m){
    ll r=0;for(;b;b>>=1,a=(a+a)%m)if(b&1)r=(r+a)%m;
    return r;
}
ll qpow(ll a,ll b,ll m){
    ll r=1;for(;b;b>>=1,a=qmul(a,a,m)%m)if(b&1)r=qmul(r,a,m)%m;
    return r;
}
int gcd(ll a,ll b){return b?gcd(b,a%b):a;}

bool check(ll x,ll u,int t,ll mo){
    ll n=qpow(x,u,mo);
    if(n==1 || n==mo-1)return 1;
    for(int i=1;i<=t;i++){
        if(n==mo-1)return 1;
        n=qmul(n,n,mo);
        if(n==1)return 0;
    }return 0;
}
bool miller_rabin(ll x){
    if(x<2)return 0;
    if(x==2)return 1;
    if(x%2==0)return 0;
    ll u=x-1,k;int s=0;
    while(u%2==0)u/=2,++s;
    for(int i=1;i<=10;i++){
        do k=R()%x;while(k==0);
        if(!check(k,u,s,x))return 0;
    }return 1;
}
ll pollard_rho(ll x,ll c){
    int i=1,k=2;
    ll x0=R()%x,y=x0;
    for(;;){
        ++i;
        x0=(qmul(x0,x0,x)+c)%x;
        ll d=gcd(abs(x0-y),x);
        if(d!=1 && d!=x)return d;
        if(x0==y)return x;
        if(i==k)k<<=1,y=x0;
    }
}
void sep(ll x){
    if(x==1)return;
    if(miller_rabin(x)){
        p[++t]=x;return;
    }
    ll d=pollard_rho(x,R()%x+1);
    sep(d);sep(x/d);
}

void dg(int w,ll d){
    if(w>t){
        if(ans>d && qpow(n,d,lcm)==1%lcm)ans=d;
        return;
    }ll m=1;
    for(int i=0;i<=ap[w];++i,m*=p[w])dg(w+1,d*m);
}

int main(){
    freopen("pow.in","r",stdin);
    freopen("pow.out","w",stdout);
    srand((int)time(0));
    scanf("%lld",&n);ans=n;
    sep(n);sort(p+1,p+t+1);
    for(int i=1;i<t;i++)if(p[i]==p[i+1]){
        printf("-1");return 0;
    }lcm=1;
    for(int i=1;i<=t;i++)lcm=lcm*(p[i]-1)/gcd(lcm,p[i]-1);
    if(gcd(lcm,n)>1){
        printf("-1");return 0;
    }t=0;sep(lcm);sort(p+1,p+t+1);ll phi=1;
    for(int i=1;i<=t;i++)if(p[i]==p[i-1])phi*=p[i];else phi*=p[i]-1;
    t=0;sep(phi);sort(p+1,p+t+1);
    for(int i=1,x=t;i<=x;i++)if(i==1)ap[t=1]=1;else if(p[i]==p[i-1])++ap[t];else ap[++t]=1,p[t]=p[i];
    dg(1,1);printf("%lld",ans);
    fclose(stdin);fclose(stdout);
    return 0;
}