BSGS算法及拓展

https://www.zybuluo.com/ysner/note/1299836

定义

一种用来求解高次同余方程的算法。
一般问题形式：求使得\(y^x\equiv z(mod\ p)\)的最小非负\(x\)。

\(BSGS\)算法

要求\(p\)是质数。
由费马小定理可知，\(y^{p-1}\equiv1(mod\ p)\)，所以暴力枚举只要枚举到\(p−1\)即可。
但是由于\(p\)一般都很大，所以一般都跑不动。。。

优化算法\(ing...\)
现在令\(x=mi−j\)(其中\(m=\lceil\sqrt p\rceil\)）。
则方程可化为\(y^{mi-j}\equiv z(mod\ p)\),
\(y^{mi}\equiv y^jz(mod\ p)\)

然后可以发现\(j<m\)（否则\(x\)就是负数）
所以我们可以暴力枚举\(j\),与所得\(y^jz(mod\ p)\)的存在哈希表里，然后再暴力枚举\(i\),最后得出结果。

还要注意一些边界：

• \(y!=0\)
• \(z=1\)时\(puts("no\ solution")\)
• \(i\)的边界是\([1,m+1]\)

一道Poj上的板子题
[SDOI2011]计算器

struct Hash_Table
{
  int h[N],cnt;
  struct Edge{int u,v,nxt;}e[N*10];
  il void clear(){memset(h,-1,sizeof(h));cnt=0;}
  il void add(re int u,re int v,re int w){e[++cnt]=(Edge){w,v,h[u]};h[u]=cnt;}
  il int Query(re int x)
  {
    re int t=x%mod;
    for(re int i=h[t];i+1;i=e[i].nxt)
      if(e[i].u==x) return e[i].v;
    return -1;
  }
  il void solve(re int y,re int z,re int p)
  {
    y%=p;z%=p;
    if(!y) {puts("no solution");return;}
    if(z==1) {puts("0");return;}
    re int m=sqrt(p)+1;clear();
    for(re int i=0,t=z;i<m;++i,t=1ll*t*y%p) add(t%mod,i,t);
    for(re int i=1,tt=ksm(y,m,p),t=tt;i<=m+1;++i,t=1ll*t*tt%p)
      {
	re int j=Query(t);if(j==-1) continue;
	printf("%d\n",i*m-j);return;
      }
    puts("no solution");
  }
}BSGS;
int main()
{
  re int y,p,z;
  while(scanf("%d%d%d",&p,&y,&z)!=EOF)
    {
      BSGS.solve(y,z,p);
    }
  return 0;
}

拓展\(BSGS\)算法

不要求\(p\)是质数。
那就说明很可能\(gcd(y,p)!=1\)，不满足费马小定理。
费马小定理提供了枚举上限，没有它这种问题就不好做了。。。

想想怎么把\(y,p\)约分。
令\(t=gcd(y,p)\)。
把方程改写成等式形式：$$y^x+kp=z$$
分析一下，可以发现\(z\)一定是\(t\)的倍数。
除\(t\):$$\frac{y}{t}y^{x-1}+\frac{p}{t}k=\frac{z}{t}$$
接下来再次检查\(gcd(y,\frac{z}{t})\)是否为\(1\)，若否，说明还可以继续约分，理由同上。

最后形式为（那个\(t\)反正是个正整数）$$\frac{y^k}{t}y\equiv\frac{z}{t}(mod\ \frac{p}{t})$$

注意边界：

• 如果\(t>1\)并且\(z\%t>0\)，方程无解
• 约分完的石子带到普通\(BSGS\)中时要带系数

咕谷模板题

#include<iostream>
#include<cmath>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<algorithm>
#define il inline
#define re register
#define ll long long
#define fp(i,a,b) for(re int i=a;i<=b;++i)
#define fq(i,a,b) for(re int i=a;i>=b;--i)
using namespace std;
const int N=5e4,mod=45807;
il ll ksm(re ll S,re ll n,re int p)
{
  re ll T=S;S=1;
  while(n)
    {
      if(n&1) S=S*T%p;
      T=T*T%p;
      n>>=1;
    }
  return S;
}
struct Hash_Table
{
  int h[N],cnt;
  struct Edge{int u,v,nxt;}e[N*10];
  il void clear(){memset(h,-1,sizeof(h));cnt=0;}
  il void add(re int u,re int v,re int w){e[++cnt]=(Edge){w,v,h[u]};h[u]=cnt;}
  il int Query(re int x)
  {
    re int t=x%mod;
    for(re int i=h[t];i+1;i=e[i].nxt)
      if(e[i].u==x) return e[i].v;
    return -1;
  }
  il void solve(re int y,re int z,re int p)
  {
    if(z==1) {puts("0");return;}
    re int k=0,a=1;
    while(233)
      {
    re int t=__gcd(y,p);if(t==1) break;
    if(z%t) {puts("No Solution");return;}
    z/=t;p/=t;++k;a=1ll*a*y/t%p;
    if(z==a) {printf("%d\n",k);return;}//有意思的地方
      }
    re int m=sqrt(p)+1;clear();
    for(re int i=0,t=z;i<m;++i,t=1ll*t*y%p) add(t%mod,i,t);
    for(re int i=1,tt=ksm(y,m,p),t=1ll*a*tt%p;i<=m+1;++i,t=1ll*t*tt%p)
      {
    re int j=Query(t);if(j==-1) continue;
    printf("%d\n",i*m-j+k);return;
      }
    puts("No Solution");
  }
}BSGS;
int main()
{
  re int y,p,z;
  while(scanf("%d%d%d",&y,&p,&z))
    {
      if(!p&&!y&&!z) break;
      BSGS.solve(y,z,p);
    }
  return 0;
}