参考资料

原根

离散对数

求质数的原根

OI-WIKI 

任意模数的BSGS算法证明

拓展欧几里得求通解

BSGS变形

 

原根

如果g是m的原根,对于任意一个数x(x<m),都可以找到一个I(x) 小于等于 φ(m),使得 gI(x) = x ,I(x)称为x的指标。也就是说在整数模m乘法群(所有与 m互素的正整数构成的等价类构成的乘法群)中,g是这个乘法群的生成元。它的生成元的个数就是它的可逆元个数,即φ(φ(m))个。也就是说,m有φ(φ(m))个原根

 除了直接运算以外,至今还没有一个办法可以找到模特定m时的原根,但假如已知模m有一个原根,则可找出它其他的原根。

离散对数

从表面通俗地理解,和普通对数一样,两边同时取log,这样就可以将幂次变成多项式,在离散对数中,做了这样的变换后,模m则变成了模φ(m)。于是我们就可以将一些乘法操作变成加法操作,来简化或者解决一些运算问题.

求质数p的最小原根

step1 : 筛素数

step2 : 分解φ(p)

step3 : 按定义判断φ(p)每个质因数d是不是原根

bool not_prime[N+5];
int prime[N+5],tot_prime,pr[N+5],tot_pr;

ll ksm(ll a,ll b,ll mod){
    ll res=1;
    while(b){
    if(b&1) res = (res*a)%mod;
    a = (a*a)%mod;
    b>>=1;
    }
    return res;
}

void sss(){ // 线性筛
    for(int i=2;i<N;i++){
    if(!not_prime[i]) prime[++tot_prime]=i;
    for(int j=1;j<=tot_prime;j++){
        int t = i*prime[j];if(t>=N) break;not_prime[t]=true;
        if(i%prime[j]==0)break;
    }
    }
}
void get(int p){ // 分解phi(p)的质因数存到pr[]中
    int temp=p-1;
    for(int i=1;i<=tot_prime&&prime[i]*prime[i]<=temp;i++)
    if(temp%prime[i]==0){
        while(temp%prime[i]==0) temp /= prime[i];
        pr[++tot_pr] = prime[i];
    }
    if(temp>1) pr[++tot_pr]=temp;
}
bool check(int d,int p){ //判断d是不是p的原根
    if (std::__gcd(d,p)!=1)return 0;
    for(int i=1;i<=tot_pr;i++)if(ksm(d,(p-1)/pr[i],p)==1)return 0;
    return 1; 
}
int get_g(int p) {// 求p的原根,p是质数
    int g=1;
    sss();get(p);
    for(int i=2;i<p;i++) if(check(i,p)) {g=i;break;}
    return g;
}

 

 

BSGS算法

给定a,b,p, 求解 ,p是个质数

时间复杂度 O(√p)

由上面的原根概念可知,这里需要求b的指标I(b)

const int MOD=76543;
ll hs[MOD],head[MOD],next[MOD],id[MOD],top;
void insert(ll x,int y){
    ll k=x%MOD;
    hs[++top]=x,id[top]=y,next[top]=head[k],head[k]=top;
}
ll find(int x){
    ll k=x%MOD;
    for(int i=head[k];i!=0;i=next[i])
    if(hs[i]==x)return id[i];
    return -1;
}
ll bsgs(ll a,ll n,ll b) {
    if(b==1)return 0;
    int m=sqrt(1.0*n),j;
    ll x=1,p=1,i;
    for(i=0;i<m;i++,p=p*a%n)insert(p*b%n,i);
    for(i=m;;i+=m){
    if((j=find(x=x*p%n))!=-1) return i-j;
    if(i>n) break;
    }
    return -1;
}

 

 

扩展BSGS算法

给定a,b,p, 求解 但是不保证p是质数

#include<bits/stdc++.h>
typedef long long ll;
const int N=1e5;
ll p,k,a;
bool not_prime[N+5];
int prime[N+5],tot_prime,pr[N+5],tot_pr;
const int MOD=76543;
ll hs[MOD],head[MOD],next[MOD],id[MOD],top;
ll ksm(ll a,ll b,ll mod){
    ll res=1;
    while(b){
    if(b&1) res = (res*a)%mod;
    a = (a*a)%mod;
    b>>=1;
    }
    return res;
}
void insert(ll x,int y){
    ll k=x%MOD;
    hs[++top]=x,id[top]=y,next[top]=head[k],head[k]=top;
}
ll find(int x){
    ll k=x%MOD;
    for(int i=head[k];i!=0;i=next[i])
    if(hs[i]==x)return id[i];
    return -1;
}
ll bsgs(ll a,ll n,ll b) {
    if(b==1)return 0;
    int m=sqrt(1.0*n),j;
    ll x=1,p=1,i;
    for(i=0;i<m;i++,p=p*a%n)insert(p*b%n,i);
    for(i=m;;i+=m){
    if((j=find(x=x*p%n))!=-1) return i-j;
    if(i>n) break;
    }
    return -1;
}
void sss(){ // 线性筛
    for(int i=2;i<N;i++){
    if(!not_prime[i]) prime[++tot_prime]=i;
    for(int j=1;j<=tot_prime;j++){
        int t = i*prime[j];if(t>=N) break;not_prime[t]=true;
        if(i%prime[j]==0)break;
    }
    }
}
void get(){ // 分解phi(p)的质因数存到pr[]中
    int temp=p-1;
    for(int i=1;i<=tot_prime&&prime[i]*prime[i]<=temp;i++)
    if(temp%prime[i]==0){
        while(temp%prime[i]==0) temp /= prime[i];
        pr[++tot_pr] = prime[i];
    }
    if(temp>1) pr[++tot_pr]=temp;
}
bool check(ll d){
    if (std::__gcd(d,p)!=1)return 0;
    for(int i=1;i<=tot_pr;i++)if(ksm(d,(p-1)/pr[i],p)==1)return 0;
    return 1; 
}
ll exgcd(ll m,ll &x,ll n,ll &y){
    ll x1=0,y1=1,x0=1,y0=0;
    ll r=(m%n+n)%n;ll q=(m-r)/n;
    x=0;y=1;
    while(r){
    x=x0-q*x1;y=y0-q*y1;
    x0=x1;y0=y1;x1=x;y1=y;
    m=n;n=r;r=m%n;
    q=(m-r)/n;
    }
    return n;
}
void solve(ll p,ll k,ll a){
    std::set<int>ans;
    int g=1;
    sss();get();
    for(int i=2;i<p;i++)if(check(i)){g=i;break;}
    ll Ia=bsgs(g,p,a);
    ll ar,br;
    ll gcd=exgcd(k,ar,p-1,br);
    if(Ia%gcd==0){
        ll sp=(p-1)/gcd;
        ar=((ar*Ia/gcd%sp)%sp+sp)%sp;
        if(ar==0)ar=sp;
        for(int i=ar;i<p;i+=sp)ans.insert(ksm(g,i,p));
    }
    int m = ans.size();
    printf("%d\n",m);
    for(std::set<int>::iterator it=ans.begin();it!=ans.end();it++)
    printf("%d\n",(*it));
}
int main(){
    scanf("%lld%lld%lld",&p,&k,&a);
    solve(p,k,a);
}

 

 BSGS算法进阶

求解xk ≡ a(mod p)   p是一个质数

该模型可以通过一系列的转化为成BSGS算法的常规模型

gI(a) ≡ a (mod p) BSGS解得I(a)

gI(x) ≡ x (mod p)

原方程变形为 (gI(x))k ≡ gI(a) (mod p)

两边取离散对数 得到 k* I(x) ≡ I(a) (mod φ(p))

即 k*I(x) + φ(p) * yy = I(a) 用扩展欧几里得算法解得一个I(x)

如何获得所有解呢?

对于关于x,y的方程a * x + b * y = gcd 来说,让x增加b/gcd,让y减少a/gcd,等式两边还相等。

对于方程a * x + b * y = c 来说,如果c % gcd(a,b) != 0,则无解。
如果c % gcd(a,b) == 0 且 c / gcd(a,b) = t,那么求出方程 a * x + b * y = gcd(a,b)的所有解x,y,将x,y乘上t,对应的x’,y’即是方程a * x + b * y = c的解

每得到一个I(x) 计算出x ,I(x)≤φ(p)

#include<bits/stdc++.h>
typedef long long ll;
const int N=1e5;
ll p,k,a;
bool not_prime[N+5];
int prime[N+5],tot_prime,pr[N+5],tot_pr;
const int MOD=76543;
ll hs[MOD],head[MOD],next[MOD],id[MOD],top;
ll ksm(ll a,ll b,ll mod){
    ll res=1;
    while(b){
    if(b&1) res = (res*a)%mod;
    a = (a*a)%mod;
    b>>=1;
    }
    return res;
}
void insert(ll x,int y){
    ll k=x%MOD;
    hs[++top]=x,id[top]=y,next[top]=head[k],head[k]=top;
}
ll find(int x){
    ll k=x%MOD;
    for(int i=head[k];i!=0;i=next[i])
    if(hs[i]==x)return id[i];
    return -1;
}
ll bsgs(ll a,ll n,ll b) {
    if(b==1)return 0;
    int m=sqrt(1.0*n),j;
    ll x=1,p=1,i;
    for(i=0;i<m;i++,p=p*a%n)insert(p*b%n,i);
    for(i=m;;i+=m){
    if((j=find(x=x*p%n))!=-1) return i-j;
    if(i>n) break;
    }
    return -1;
}
void sss(){ // 线性筛
    for(int i=2;i<N;i++){
    if(!not_prime[i]) prime[++tot_prime]=i;
    for(int j=1;j<=tot_prime;j++){
        int t = i*prime[j];if(t>=N) break;not_prime[t]=true;
        if(i%prime[j]==0)break;
    }
    }
}
void get(){ // 分解phi(p)的质因数存到pr[]中
    int temp=p-1;
    for(int i=1;i<=tot_prime&&prime[i]*prime[i]<=temp;i++)
    if(temp%prime[i]==0){
        while(temp%prime[i]==0) temp /= prime[i];
        pr[++tot_pr] = prime[i];
    }
    if(temp>1) pr[++tot_pr]=temp;
}
bool check(ll d){
    if (std::__gcd(d,p)!=1)return 0;
    for(int i=1;i<=tot_pr;i++)if(ksm(d,(p-1)/pr[i],p)==1)return 0;
    return 1; 
}
ll exgcd(ll m,ll &x,ll n,ll &y){
    ll x1,y1,x0,y0;
    x0=1;y0=0;
    x1=0;y1=1;
    ll r=(m%n+n)%n;
    ll q=(m-r)/n;
    x=0;y=1;
    while(r){
    x=x0-q*x1;
    y=y0-q*y1;
    x0=x1;y0=y1;
    x1=x;y1=y;
    m=n;n=r;r=m%n;
    q=(m-r)/n;
    }
    return n;
}
void solve(ll p,ll k,ll a){
    std::set<int>ans;
    int g=1;
    sss();get();
    for(int i=2;i<p;i++)if(check(i)){g=i;break;}
    ll Ia=bsgs(g,p,a);
    ll ar,br;
    ll gcd=exgcd(k,ar,p-1,br);
    if(Ia%gcd==0){
        ll sp=(p-1)/gcd;
        ar=((ar*Ia/gcd%sp)%sp+sp)%sp;
        if(ar==0)ar=sp;
        for(int i=ar;i<p;i+=sp)ans.insert(ksm(g,i,p));
    }
    int m = ans.size();
    printf("%d\n",m);
    for(std::set<int>::iterator it=ans.begin();it!=ans.end();it++)
    printf("%d\n",(*it));
}
int main(){
    scanf("%lld%lld%lld",&p,&k,&a);
    solve(p,k,a);
}