【JZOJ5745】幂

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

Input
Inp

Output
Out

Sample

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

Data Constraint
Data


analysis

anka(modn)∀ank≡a(modn)

很容易发现当n有非1的平方因子时一定无解,因为设n=p2q(p>1)n=p2q(p>1),当a=pa=p时,又因为n2n≥2,那么肯定找不到k
那么有n=pin=∏pi
该式子可以分解为多个anka(modpi)ank≡a(modpi)
则会有nk1(modpi1)nk≡1(modpi−1)
因此,设L=LCM(p11,p21,,ps1)L=LCM(p1−1,p2−1,⋯,ps−1)
那么
nk1(modL)nk≡1(modL)

同时当gcd(n,L)>1gcd(n,L)>1也无解

nφ(L)1(modL)∵nφ(L)≡1(modL)

k|φ(L)k|φ(L)

枚举φ(L)φ(L)的因数即可
可能用到快速幂和快速加。

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

如何快速质因数分解呢?

这里需要用到Miller_RabinMiller_RabinPollard_rhoPollard_rho

Miller_RabinMiller_Rabin快速判定质数法
说到质数,很容易想到费马小定理ap11(modp)     (1a<p)ap−1≡1(modp)     (1≤a<p)
难道用这个判就好了吗?
很可惜不能,有一些强伪素数对于任意的a都满足这个式子,但是他是个合数,如561
但是有一个判定质数的充分必要条件,那就是1<a<p1a2≢1(modp)∀1<a<p−1a2≢1(modp)
那么,令p1=2uvp−1=2uv
只用随机a,判断ap11(modp)ap−1≡1(modp)
a2uv1(modp)a2u′v≡1(modp)同时a2uv±1(modp)a2u′v≡±1(modp)
一次随机的正确率为3434则随机个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_rhoPollard_rho快速找因数
首先先引入生日悖论
若要从1~1000中选取1,求选取成功的几率?
答案是1100011000
但是如果变为随机两个数,求其差为1的概率为多少,这就变成了15001500
如果取三个数,存在两者的差为1的概率,将增大为35003500
通过实验,若选30个数概率将超过50%50%
重复实验证明,当1~n中选择k个数,使得存在两个数差值为某值,若knk≤n则其概率将大于50%50%

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

Pollard’s rho算法
其实并不需要将k个数弄出来,只需一个一个地生成并检查连续的两个数,反复执行这个步骤并希望能够得到我们想要的数。
能够符合这个条件的变换则是伪随机f(x)=(x2+c)modnf(x)=(x2+c)modn
也就是xi+1=f(xi)xi+1=f(xi)fxfx和另外一个x进行判断,而另外一个也要求能通过之前的变换而来
这种方案最终是存在循环的,如何判断循环呢
要个很巧妙的方法
这里拿x2i+12i+11x2i+1→2i+1−1x2i+1x2i+1进行判断,若判断的过程中,存在x2i+k=x2i+1x2i+k=x2i+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;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值