[Sgu 261]Discrete Roots

Time Limit : 1s Memory Limit : 64MB

题目大意

给出三个整数p,k,a,其中p为质数，求出所有满足x^k=a (mod p),0<=x<=p-1的x。

解题报告

首先假设我们都学会了BSGS， 原根， 扩展欧几里得。。。 原根看这里 看到这一题最先想到BSGS， 毕竟公式很像，但是求得是底数而不是指数。

看题目给出的$$x^k \equiv a \pmod p $$ 暴力枚举一发求出原根r，那么有唯一的一个i满足$$r^i \equiv x \pmod p$$ 也有唯一的一个j满足$$r^j \equiv a \pmod p$$ 由于已知r,a,p,这就是一个BSGS的模型了，求离散对数即可 然后现在要求i就求出了x 然后从题目原式结合一下得到$$r^{ik} \equiv a \pmod p$$ 把a换成r^j(mod p下)得到$$r^{ik} \equiv r^j \pmod p$$ 因为r是p的原根， 那么一定有$$ i k \equiv j \pmod {p-1}$$ 这是个明显的线性同余方程， 用扩展欧几里得求解就可以啦^.^

需要注意的是p-1不一定和k互质，需要除一下。

还有一点， 如果在bzoj上提交此题，每个答案之间要换行，不然会PE = =！ 接下来是AC代码

代码

#include <cstdio>
#include <cstring>
#include <cstdlib>
#include <algorithm>
#include <cmath>
#include <iostream>

#include <map>

typedef long long LL;
const int MAXN = 110000;
std :: map<LL, LL> mp;

void exgcd(LL a, LL b, LL & x, LL & y, LL & d) {
    if(b == 0) {
        x = 1; y = 0; d = a;
        return ;
    }
    exgcd(b, a%b, x, y, d);
    LL t = x;
    x = y;
    y = t - a/b*y;
}
//A^x
void setup (LL A, LL C) {
    LL i = 1, u = ceil(sqrt(C)+0.5);
    for(LL r = A%C; i<u; i++, (r*=A)%=C) mp[r] = i;
    mp[1] = 0;
}

LL advPow (LL a, int b, LL c) {
    a %= c;
    LL ret = 1;
    while(b) {
        if(b & 1) ret *= a, ret %= c;
        a *= a;
        a %= c;
        b >>= 1;
    }
    return ret;
}

//A^x = B (mod C)
LL BSGS (LL A, LL B, LL C) {
    if(A % C == 0)
        return 0;
    setup(A, C);
    int r = ceil(sqrt(C+0.5)); // иох║уШ
    LL D = 1, T = advPow(A, r, C);
    for(int i = 0; i<r; i++) {
        LL a, b, t;
        exgcd(D, C, a, b, t);
        LL c = C/t;
        a = (a*B/t%c+c)%c;
        //prLLf("a:%d i:%d\n", a, i);
        if(mp.find(a) != mp.end())
                return mp[a] + i*r;
        D = D*T%C;
    }
    return -1;
}

const int PMAX = 110000;
LL prime[PMAX/10];
int ptot;
bool pvis[PMAX];

void phi () {
    int r = sqrt(PMAX+0.5);
    for(int i = 2; i<=r; i++)
        if(!pvis[i]) for(LL j = i+i; j<PMAX; j+=i) pvis[j] = true;
    for(int i = 2; i<PMAX; i++)
        if(!pvis[i]) prime[++ptot] = i;
}

//r^a = b(mod P)
LL divn[30];
int getOriginRt (LL P) { // P must be a prime, or eular(P) != P-1
    LL up = sqrt((P-1)+0.5);
    LL dtot = 0;
    LL t = P-1;
    for(int i = 1; prime[i]<up; i++)
        if(t % prime[i] == 0) {
            divn[++dtot] = prime[i];
            while(t % prime[i] == 0) t /= prime[i];
        }
    for(int r = 2; r<P; r++) {
        if(advPow(r, P-1, P) != 1) continue ;
        bool flag = true;
        for(int i = 1; i<=dtot; i++)
            if(advPow(r, (P-1)/divn[i], P) == 1)
                {flag = false; break;}
        if(flag) return r;
    }
    return *((int*)NULL);
}

LL ans[MAXN];
int atot;

using std :: cin;
using std :: cout;

int main () {
    LL k, a, P;
    //x^k = a (mod p)
    cin >> P >> k >> a;
    if(!a)
    {
        puts("1");puts("0");
        return 0;
    }
    phi();
    LL r = getOriginRt(P);
    //prLLf("r:%d\n", r);
    LL j = BSGS(r, a, P);
    //prLLf("j:%d\n", j);
    LL x, y, g;
    exgcd(k, P-1, x, y, g);
    if(j % g) {
        puts("0");
        return 0;
    }
    j/=g;
    //printf("r:%d\n", r);
    LL c = (P-1)/g;
    x = (x*j%c+c)%c;
    while(x < P-1) {
        ans[++atot] = advPow(r, x, P);
        x += c;
    }
    std :: sort(ans + 1, ans + (atot+1));
    cout << atot << '\n';
    for(int i = 1; i<=atot; i++)
        cout << ans[i] << '\n';
}

PS : 新博客上的第一篇题解，好兴奋啊！ PPS : 装弱是会遭雷劈的！: ( @Ender