模板BSGS（SDOI2011计算器） 模板EXBSGS

BSGS和EXBSGS是OI中用于解决A^xΞB(mod C)的常用算法。

1.BSGS

BSGS用于A,C互质的情况。

令m=sqrt(C)，此时x可表示为i*m+j。

式中i和j都<=sqrt(C)

原式Ax≡B(mode C) -->Ai*m * Aj≡B(mode C)

枚举Ai*m,此时Ai*m相当于系数。//O(sqrt(C))

现在我们可用exgcd/费马小定理求逆元算出Aj%C的值

通过预处理将A1~m存入map/哈希表。//O(1)//用map会多一个log

解决了。

时间复杂度O(sqrt(C))，思想类似meet_in_the_middle。

代码：

#include<map>
#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define ll long long
ll fastpow(ll x,int y,int mod)
{
    ll ret = 1ll;
    while(y)
    {
        if(y&1)ret=ret*x%mod;
        x=x*x%mod;
        y>>=1;
    }
    return ret;
}
void ORZ(){printf("Orz, I cannot find x!\n");}
ll F2(ll y,int z,int p)
{
    ll tmp = fastpow(y,p-2,p);
    tmp = tmp*z%p;
    return tmp;
}
int hed[60050],cnt=0;
struct EG
{
    int to,w,nxt;
}e[70050];
void ae(int f,int t,int w)
{
    e[++cnt].to = t;
    e[cnt].w=w;
    e[cnt].nxt = hed[f];
    hed[f] = cnt;
}
int find(int x)
{
    int now = x%60000;
    for(int j=hed[now];j;j=e[j].nxt)
        if(e[j].to==x)
            return e[j].w;
    return -1;
}
void BSGS(ll y,int z,int p)
{
    if(!y&&z){ORZ();return ;}
    cnt=0;
    memset(hed,0,sizeof(hed));cnt=0;
    ae(1,1,0);
    int m = (int)sqrt(p);
    ll now = 1;
    for(int i=1;i<=m;i++)
    {
        now = now*y%p;
        if(find(now)==-1)
        {
            ae(now%60000,now,i);
        }
    }
    ll u = 1;
    for(int i=0;i<=m+2;i++)
    {
        ll tmp = F2(u,z,p);
        ll ans = find(tmp);
        if(~ans)
        {
            printf("%lld\n",(ans+i*m)%p);
            return ;
        }
        u=u*now%p;
    }
    ORZ();
}
int T,K,y,z,p;
int main()
{
    scanf("%d%d",&T,&K);
    while(T--)
    {
        scanf("%d%d%d",&y,&z,&p);
        y%=p;
        if(K==1)printf("%lld\n",fastpow(y,z,p));
        else if(K==2)
        {
            z%=p;
            if(!y&&z){ORZ();continue;}
            printf("%lld\n",F2(y,z,p));
        }else
        {
            z%=p;
            BSGS(y,z,p);
        }
    }
    return 0;
}

2.EXBSGS

EXBSGS用于解决C不一定是质数的问题。

对于方程Ax≡B(mode C)(gcd(A,C)!=1)：

设d=gcd(A,C)。

如果B%d!=0的话，要么B==0，要么无解。

这时我们可以将三者同时/d，得到：

(A/d)*Ax-1≡ B/d (mode C/d)

此时A和C/d依然可能不互质（比如A=5,C=25）。

循环这个过程就可以啦。

我们最后得到的方程形式为：

(Ai/∏d)*Ax-i ≡ B/∏d(mode C/∏d)

其实就是D*Ax-i≡ B’(mode C’)，然后套BSGS就好了。

代码：

#include<cmath>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define ll long long
void ORZ()
{
    printf("No Solution\n");
}
void exgcd(ll a,ll b,ll &x,ll &y)
{
    if(!b)
    {
        x=1,y=0;
        return ;
    }
    exgcd(b,a%b,y,x);
    y-=a/b*x;
}
ll inv(ll a,ll b)
{
    ll x,y;
    exgcd(a,b,x,y);
    return (x%b+b)%b;
}
ll gcd(ll x,ll y)
{
    return y?gcd(y,x%y):x;
}
ll F2(ll y,ll z,ll p)
{
    y%=p,z%=p;
    ll ret = inv(y,p);
    return ret*z%p;
}
int hed[40050],cnt;
struct EG
{
    ll to,nxt,w;
}e[100050];
void ae(ll f,ll t,ll w)
{
    e[++cnt].to = t;
    e[cnt].nxt = hed[f];
    e[cnt].w = w;
    hed[f] = cnt;
}
ll find(ll x)
{
    for(int j=hed[x%40000];j;j=e[j].nxt)
        if(e[j].to==x)
            return e[j].w;
    return -1;
}
ll BSGS(ll y,ll z,ll p)
{
    memset(hed,0,sizeof(hed));cnt=0;
    ae(1,1,0);
    y%=p,z%=p;
    if(!y&&z)
    {
        ORZ();
        return -1;
    }
    ll now = 1ll,m = (ll)sqrt(p);
    for(int i=1;i<=m;i++)
    {
        now=now*y%p;
        if(find(now)==-1)
            ae(now%40000,now,i);
    }
    ll u = 1ll;
    for(int i=0;i<=m+2;i++)
    {
        ll tmp = F2(u,z,p);
        ll ans = find(tmp);
        if(~ans)
            return ans+i*m;
        u=u*now%p;
    }
    ORZ();
    return -1;
}
void EXBSGS(ll y,ll z,ll p)
{
    if(!y&&z){ORZ();return ;}
    if(z==1){printf("0\n");return ;}
    ll d = gcd(y,p),cnt=0,D=1ll;
    while(d!=1)
    {
        if(z%d){ORZ();return ;}
        cnt++;
        z/=d,p/=d;
        D=D*(y/d)%p;
        d=gcd(y,p);
        if(D==z){printf("%lld\n",cnt);return ;}
    }
    ll ans = BSGS(y,z*inv(D,p)%p,p);
    if(ans!=-1)printf("%lld\n",(ans+cnt)%p);
}
ll y,z,p;
int main()
{
    while(scanf("%lld%lld%lld",&y,&p,&z))
    {
        if(!y&&!z&&!p)break;
        EXBSGS(y,z,p);
    }
    return 0;
}