結果

問題 No.981 一般冪乗根
ユーザー ace_amuroace_amuro
提出日時 2021-01-26 17:05:12
言語 C++14
(gcc 12.3.0 + boost 1.83.0)
結果
AC  
実行時間 4,425 ms / 6,000 ms
コード長 5,715 bytes
コンパイル時間 967 ms
コンパイル使用メモリ 85,236 KB
実行使用メモリ 96,424 KB
最終ジャッジ日時 2023-09-05 20:52:06
合計ジャッジ時間 86,693 ms
ジャッジサーバーID
（参考情報） 		judge11 / judge15
このコードへのチャレンジ
（要ログイン）

テストケース

テストケース表示
入力 結果 実行時間
実行使用メモリ
testcase_00 AC 375 ms
96,424 KB
testcase_01 AC 348 ms
49,148 KB
testcase_02 AC 387 ms
49,148 KB
testcase_03 AC 328 ms
49,144 KB
testcase_04 AC 394 ms
49,092 KB
testcase_05 AC 312 ms
49,084 KB
testcase_06 AC 385 ms
49,144 KB
testcase_07 AC 390 ms
49,096 KB
testcase_08 AC 429 ms
49,160 KB
testcase_09 AC 347 ms
49,152 KB
testcase_10 AC 382 ms
49,204 KB
testcase_11 AC 371 ms
49,156 KB
testcase_12 AC 384 ms
49,100 KB
testcase_13 AC 321 ms
49,196 KB
testcase_14 AC 360 ms
49,132 KB
testcase_15 AC 392 ms
49,272 KB
testcase_16 AC 287 ms
49,156 KB
testcase_17 AC 346 ms
49,144 KB
testcase_18 AC 317 ms
49,208 KB
testcase_19 AC 381 ms
49,160 KB
testcase_20 AC 342 ms
49,116 KB
testcase_21 AC 378 ms
49,264 KB
testcase_22 AC 385 ms
49,144 KB
testcase_23 AC 383 ms
49,116 KB
testcase_24 AC 375 ms
49,136 KB
testcase_25 AC 85 ms
18,380 KB
testcase_26 AC 84 ms
18,500 KB
testcase_27 AC 477 ms
49,268 KB
testcase_28 AC 4,425 ms
49,116 KB
evil_60bit1.txt AC 460 ms
49,188 KB
evil_60bit2.txt AC 463 ms
49,148 KB
evil_60bit3.txt AC 415 ms
49,244 KB
evil_hack AC 98 ms
49,116 KB
evil_hard_random AC 150 ms
49,132 KB
evil_hard_safeprime.txt AC 127 ms
18,480 KB
evil_hard_tonelli0 AC 226 ms
49,396 KB
evil_hard_tonelli1 AC 5,903 ms
49,212 KB
evil_hard_tonelli2 AC 3,323 ms
49,240 KB
evil_hard_tonelli3 AC 223 ms
49,216 KB
evil_sefeprime1.txt AC 125 ms
18,440 KB
evil_sefeprime2.txt AC 125 ms
18,428 KB
evil_sefeprime3.txt AC 126 ms
18,408 KB
evil_tonelli1.txt TLE -
evil_tonelli2.txt TLE -
権限があれば一括ダウンロードができます

ソースコード

diff # 

// #define TEST
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<cmath>
#include<string>
#include<algorithm>
#include<vector>
#include<queue>
#include<map>
using namespace std;
typedef long long LL;
const int M=1e7+10;
bool is_p[M];
int prime[M], primelen;
int euler(int N) {
    int len = 0;
    memset(is_p, true, sizeof(is_p));
    for(int i = 2; i < N; i++) {
        if(is_p[i]) {
            prime[len] = i;
            len++;
        }
        for(int j = 0; j < len && prime[j] <= i; j++) {
            if(i * prime[j] >= N) break;
            is_p[i * prime[j]] = false;
            if(i % prime[j] == 0) break;
        }
    }
    return len;
}

LL qmul(LL a,LL b,LL mod){
    return (LL)((__int128_t)a*b%mod);
    /*
    LL ans=0;
    while(b>0)
    {
        if(b&1) ans=(ans+a)%mod;
        a=(a+a)%mod;
        b>>=1;
    }
    return ans;
    */
}

LL qpow(LL a, LL b, LL mod){
    if(b==0)return 1ll;
    LL ans=1;
    while(b>0){
        if(b&1)ans=qmul(ans,a,mod);
        a=qmul(a,a,mod);
        b>>=1;
    }
    return ans;
}



LL gcd(LL a,LL b){
    while(b){
        LL t=a%b;
        a=b;
        b=t;
    }
    return a;
}

LL exgcd(LL a,LL b,LL &x,LL &y){
    if(b == 0){
        x = 1, y = 0;
        return b;
    }
    LL d = exgcd(b, a%b, x, y);
    LL t = y;
    y = x - y * (a / b);
    x = t;
    return d;
}

LL solveLinear(LL a,LL b,LL p){
    a=(a%p+p)%p;b=(b%p+p)%p;
    LL x,y;
    LL d=exgcd(a,p,x,y);
    if(b%d!=0) return -1;
    return (qmul(b/d,x,p/d)+p/d)%(p/d);
}

LL invP(LL a,LL p){
    a=(a%p+p)%p;
    if(a==0) return -1;
    LL x,y;
    exgcd(a,p,x,y);
    return (x%p+p)%p;
}

int factorlist(LL x,LL z,LL* p,int* k){
    int cnt=0;
    LL y=gcd(z/x,x);
    LL cap=pow(z,0.25)+1;
    for(int i=0;prime[i]<=cap;i++){
        if(y%prime[i]!=0) continue;
        cnt++;
        p[cnt]=prime[i];
        while(y%prime[i]==0){
            y/=prime[i];
        }
    }
    if(y>1){
        cnt++;
        p[cnt]=y;
    }
    for(int i=1;i<=cnt;i++){
        k[i]=0;
        while(x%p[i]==0){
            x/=p[i];
            k[i]++;
        }
    }
    p[0]=x;
    return cnt;
}

const LL hashmod=1000003;
LL phash[hashmod*2];
LL r_origin[hashmod*2];
void makehash(LL Bps,LL p,LL P){
    memset(phash,-1,sizeof(phash));
    memset(r_origin,-1,sizeof(r_origin));
    LL r=1;
    LL m=sqrt(p)+1;
    phash[r%hashmod]=0;
    r_origin[r%hashmod]=r;
    for(LL i=1;i<m;i++){
        r=qmul(r,Bps,P);
        LL k=r%hashmod;
        while(phash[k]!=-1) k++;
        phash[k] = i;
        r_origin[k]=r;
    }
}

LL bsgs(LL a, LL b, LL p, LL P){//a^r*b=1 mod P
    if(a % P == 0 || b % P == 0) return -1;
    if(b % P == 1) return 0;
    if(qmul(a, b, P) == 1) return 1;
    LL m=sqrt(p)+1;
    LL t=b;
    LL am=qpow(a,m,P);
    for(LL i=1;i<=m;i++){
        t=qmul(t,am,P);
        LL k=t%hashmod;
        if(phash[k]!=-1){
            while(phash[k]!=-1){
                if(r_origin[k]==t) return m * i - phash[k];
                k++;
            }
        }
    }
    return -1;
}

LL tonelli(LL P,LL p,int k,LL A,LL invA){
    LL n=1;
    for(int i=0;i<k;i++) n*=p;
    if(qpow(A,(P-1)/n,P)!=1) return -1;
    LL v,d,q=P-1;
    int s=0;
    while(q%p==0){q/=p;s++;}
    exgcd(n,q,v,d);
    v=(v%q+q)%q;
#ifdef TEST
    printf("n=%lld q=%lld p=%lld s=%d k=%d v=%lld\n",n,q,p,s,k,v);
#endif
    LL X=qpow(A,v,P);
    int t=k+1;
    LL B=2;
    while(qpow(B,(P-1)/p,P)==1) B++;
#ifdef TEST
    printf("W=%lld\n",B);
#endif
    B=qpow(B,q,P);
    LL invB=invP(B,P);
#ifdef TEST
    printf("B=%lld B^-1=%lld\n",B,invB);
#endif
    LL E=X;
    for(int i=1;i<=k;i++) E=qpow(E,p,P);
    E=qmul(E,invA,P);
#ifdef TEST
    printf("X=%lld E=%lld\n",X,E);
#endif
    LL Bpow[100];Bpow[0]=B;
    for(int i=1;i<=s-1;i++) Bpow[i]=qpow(Bpow[i-1],p,P);
    LL Bps=Bpow[s-1];
    makehash(Bps,p,P);
    while(E!=1){
        LL Eps=E;
        for(int i=1;i<=s-t;i++) Eps=qpow(Eps,p,P);
#ifdef TEST
        printf("t=%d Bps=%lld Eps=%lld\n",t,Bps,Eps);
#endif
        LL r=bsgs(Bps,Eps,p,P);
#ifdef TEST
        printf("r=%lld\n",r);
#endif
        E=qmul(E,qpow(Bpow[t-1],r,P),P);
        X=qmul(X,qpow(Bpow[t-k-1],r,P),P);
#ifdef TEST
        printf("X=%lld E=%lld\n",X,E);
#endif
        t++;
    }
    return X;
}


LL work(LL P,LL n,LL A){
    A=(A%P+P)%P;
    LL g=gcd(n,P-1);
    if(qpow(A,(P-1)/g,P)!=1) return -1;
    LL inv=invP(n/g,(P-1)/g);
    A=qpow(A,inv,P);
    LL p[1000],X1,X2,n1,n2;int k[1000];
    LL invA=invP(A,P);
#ifdef TEST
    printf("P=%lld n=%lld g=%lld A=%lld invA=%lld\n",P,n,g,A,invA);
#endif
    int cnt=factorlist(g,P-1,p,k);
#ifdef TEST    
    printf("p:");for(int i=0;i<=cnt;i++)printf(" %lld",p[i]);printf("\n");
    printf("k:");for(int i=0;i<=cnt;i++)printf(" %d",k[i]);printf("\n");
#endif
    n1=p[0];
    LL v=invP(n1,(P-1)/n1);
    X1=qpow(A,v,P);
    for(int i=1;i<=cnt;i++){
#ifdef TEST
        printf("------------------------\n");
#endif
        X2=tonelli(P,p[i],k[i],A,invA);
        if(X2==-1) return -1;
        n2=1;
        for(int j=1;j<=k[i];j++) n2*=p[i];
#ifdef TEST
        printf("n1=%lld X1=%lld n2=%lld X2=%lld\n",n1,X1,n2,X2);
#endif
        LL alpha,beta;
        exgcd(n1,n2,alpha,beta);
        if(alpha<=0){
            X1=qpow(X1,beta,P);
            X2=qpow(invP(X2,P),-alpha,P);
        }
        else if(beta<=0){
            X1=qpow(invP(X1,P),-beta,P);
            X2=qpow(X2,alpha,P);      
        }
        X1=qmul(X1,X2,P);
        n1*=n2;
    }
    return X1;

}

int main(){
    primelen=euler(M);
    int T;
    scanf("%d",&T);
    for(;T>0;T--){
        LL P,n,A,X;
        scanf("%lld%lld%lld",&P,&n,&A);
        X=work(P,n,A);
        printf("%lld\n",X);
    }
    return 0;
}
0