[省選前題目整理][SGU 261]Discrete Roots(擴充歐幾裡得+中國剩餘定理+原根+大步小步算法)

題目連結

http://acm.hust.edu.cn/vjudge/problem/viewProblem.action?id=22207

題目大意

已知 k、y、p ，求 xk≡y(modp) 的所有可行解

思路

首先，我們求出 modp 的原根 g 。

不妨設gq≡y(modp)， q 可以通過大步小步算法求得。

由于在1≤x<p時， modp 意義下所有的數字均可以用 g 的幂次表示。是以一定存在等式gi≡k(modp−1)。這個可以通過中國剩餘定理找出所有的可行解。

—————–以下來閑扯一下大步小步算法——————-

大步小步算法看上去非常的暴力，個人認為其思路就是基于分塊思想。假設已知 x、y、p ，求 xk≡y(modp)的k ，那麼首先我們要對于所有的 i≤p√ ,在map中記錄下 xi 對應于 i 的映射。然後維護一個值xt，初始時 t=0 ，每次 t+=p√ ,若在map中能夠查找到 yxt(m−2) (這個操作複雜度隻是對數),那麼就找到了答案， k=t+mp[yxt(m−2)] ， mp[yxt(m−2)]=yxt(m−2)在 map中映射的幂

代碼

#include <iostream>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <algorithm>
#include <vector>
#include <map>
#include <cmath>

using namespace std;

typedef long long int LL;

vector<LL>divisor,ans;

LL extGCD(LL a,LL b,LL &x,LL &y)
{
    if(!b)
    {
        x=;
        y=;
        return a;
    }
    LL gcd=extGCD(b,a%b,x,y);
    LL tmp=x;
    x=y;
    y=tmp-a/b*y;
    return gcd;
}

LL fastPow(LL base,LL pow,LL mod)
{
    LL ans=;
    while(pow)
    {
        if(pow&) ans=(ans*base)%mod;
        base=(base*base)%mod;
        pow>>=;
    }
    return ans;
}

bool g_test(LL g,LL p) //檢查g是否是mod p的原根
{
    for(int i=;i<(int)divisor.size();i++)
        if(fastPow(g,(p-)/divisor[i],p)==)
            return false;
    return true;
}

LL PrimitiveRoot(LL p) //求mod p的原根
{
    divisor.clear();
    LL tmp=p-;
    for(LL i=;i*i<=tmp;i++)
    {
        if(tmp%i==)
        {
            divisor.push_back(i);
            while(tmp%i==) tmp/=i;
        }
    }
    if(tmp!=) divisor.push_back(tmp);
    LL g=; //爆枚原根g
    while(++g)
        if(g_test(g,p))
            return g;
    return -;
}

LL DiscreteLog(LL x,LL n,LL m) //求離散對數,x^y=n(mod m)，求y
{
    map<LL,LL>table; //一個有序表
    LL sqrtm=(LL)(sqrt(m)+);
    LL tmp=; //x^i=tmp
    for(int i=;i<sqrtm;i++)
    {
        table[tmp]=i;
        tmp=tmp*x%m;
    }
    LL step=tmp; //t=x^sqrt(m)
    tmp=;
    for(int i=;i<sqrtm;i++)
    {
        LL d=n*fastPow(tmp,m-,m)%m;
        if(table.count(d))
            return i*sqrtm+table[d];
        tmp=tmp*step%m;
    }
    return -;
}

void CRT(LL a,LL b,LL n) //ax=b(mod n)
{
    LL x,y,gcd;
    gcd=extGCD(a,n,x,y);
    if(b%gcd==) //有解
    {
        x=(x%n+n)%n;
        ans.push_back(x*(b/gcd)%(n/gcd));
        for(LL i=;i<gcd;i++)
            ans.push_back((ans[]+i*n/gcd)%n);
    }
}

int main()
{
    LL a,p,g,q,k;
    while(scanf("%lld%lld%lld",&p,&k,&a)!=EOF) //x^k=a(mod p)
    {
        if(!a) //a=0要特判,x隻能為0
        {
            printf("1\n0\n");
            continue;
        }
        g=PrimitiveRoot(p); //g=mod p的原根
        q=DiscreteLog(g,a,p); //g^q=a(mod p)
        CRT(k,q,p-); //問題變為kx=q(mod p-1)的所有正整數解
        for(int i=;i<(int)ans.size();i++)
            ans[i]=fastPow(g,ans[i],p); //ans[i]=g^x mod p
        sort(ans.begin(),ans.end());
        printf("%d\n",(int)ans.size());
        for(int i=;i<(int)ans.size();i++) printf("%lld\n",ans[i]);
    }
    return ;
}