BSGS/exBSGS 学习笔记【同余】
点击量：317
〇、前言
BSGS 是用来求形如 $a^x\equiv b\pmod p$ 的最小非负整数解的算法。
当 $\gcd(a,p)=1$ 时可以直接使用 BSGS，当 $\gcd(a,p)\ne 1$ 时需要用到 exBSGS 进行转化。
参考资料：
PoPoQQQ 《原根与指标》
ButterflyDew 数论ex
一、分块
对于 $a^x$，在 $\text{mod}$ 质数 $p$ 的意义下有 $\{1,2,\ldots,p-1\}$ 共 $p-1$ 种取值。
考虑把 $a^0,a^1,\ldots,a^{p-2}$ 分成 $\left\lceil\frac p{\lfloor\sqrt p\rfloor}\right\rceil$ 段，每一段的长度为 $\lfloor\sqrt p\rfloor$。预处理出 $a^1,a^2,\ldots,a^{\lfloor\sqrt p\rfloor}$ 和 $a^{\sqrt p},a^{2\sqrt p},\ldots,a^{\left\lceil\frac p{\lfloor\sqrt p\rfloor}\right\rceil \sqrt p}$，就可以用这些元素组成 $a^0,a^1,\ldots,a^{p-2}$ 中的任意一个数了。
即把 $x$ 表示成 $x=k\lfloor\sqrt p\rfloor-r$ 的形式，减法可以防止后面求逆元的麻烦。
二、meet-in-middle
在 $a^x\equiv b\pmod p$ 中，我们已知 $a,b,p$，因此 $a^x$ 的值是已知的。
在上面一段中，把 $x$ 表示成 $x=k\lfloor\sqrt p\rfloor-r$，则 $a^{k\lfloor\sqrt p\rfloor-r}\equiv b\pmod p$。
移项得 $a^{k\lfloor\sqrt p\rfloor}\equiv a^rb\pmod p$。
如果我们把 $a^1b,a^2b,\ldots,a^{\lfloor\sqrt p\rfloor}b$ 放在哈希表 $A$ 中，然后遍历 $a^{\lfloor\sqrt p\rfloor},a^{2\lfloor\sqrt p\rfloor},\ldots,a^{\left\lceil\frac p{\lfloor\sqrt p\rfloor}\right\rceil \sqrt p}$。
可以看出，等式右边已经被存在表中了，现在考虑等式左边。如果存在 $a^x\equiv b\pmod p$，那么就一定存在一个 $k$，可以在哈希表中找到和 $a^{k\lfloor\sqrt p\rfloor}$ 相等的元素。
现在要解决的问题是找到最小非负整数解。
对于大块来说，$a^{\lfloor\sqrt p\rfloor},a^{2\lfloor\sqrt p\rfloor},\ldots,a^{\left\lceil\frac p{\lfloor\sqrt p\rfloor}\right\rceil \sqrt p}$ 从小到大枚举就可以了，在哈希表中查到的顺序就不一定了。
事实上，我们用链表或前向星作为哈希表，从小到大插入 $a^1,a^2,\ldots,a^{\lfloor\sqrt p\rfloor}$。在同一链表中如果有多个符合条件的，那么就选择最后面一个；在同一前向星表中如果有多个符合条件的，那么取最前面一个。
三、BSGS
当 $\gcd(a,p)=1$ 时，$a^0,a^1,\ldots,a^{p-2}$ 刚好是 $p-1$ 个不同的数，此时直接用上面的 meet-in-middle 就可以解决了。
注意一点，我的代码风格习惯把在小段把 $a$ 滚动到 $a^{\lfloor\sqrt p\rfloor}$ 然后直接用。此时如果滚动的是 $a^ib$，那么后面还要除以 $b$。这时 $b$ 可能没有模 $p$ 的逆元（如 $b=0$），因此只滚动 $a^i$ 比较好，每次手动乘 $b$ 不会很慢。
注意一点，我的代码风格习惯把在小段把 $a$ 滚动到 $a^{\lfloor\sqrt p\rfloor}$ 然后直接用。此时如果滚动的是 $a^ib$，那么后面还要除以 $b$。这时 $b$ 可能没有模 $p$ 的逆元（如 $b=0$），因此只滚动 $a^i$ 比较好，每次手动乘 $b$ 不会很慢。
代码 of POJ 2417
#include<cstdio>
#include<cstring>
#include<map>
#include<cstdlib>
#include<cmath>
#define ll long long
#define gc() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
char buf[2100000],*p1=buf,*p2=buf;
int read()
{
int x=0;
char ch=gc();
while(ch<'0'||ch>'9')
{
if(ch==EOF)
exit(0);
ch=gc();
}
while(ch>='0'&&ch<='9')
{
x=x*10+(ch&15);
ch=gc();
}
return x;
}
ll qpow(ll x,ll y,ll p)
{
ll ans=1;
while(y)
{
if(y&1)
ans=ans*x%p;
x=x*x%p;
y>>=1;
}
return ans;
}
struct edge//前向星
{
int n,nxt,x,v;
}e[100000];
int head[10007],ecnt=-1,cnt=0;
void add(int bs,int x,int v)
{
e[++ecnt]=(edge){++cnt,head[bs],x,v};
head[bs]=ecnt;
}
int Find(int n,int x)//找到符合条件的第一个就可以返回答案了
{
for(int i=head[n];~i;i=e[i].nxt)
if(e[i].x==x)
return e[i].v;
return 0;
}
int main()
{
#ifdef wjyyy
freopen("a.in","r",stdin);
#endif
int p,b,n;
while(1)
{
for(int i=0;i<10007;++i)
head[i]=-1;
ecnt=-1;
p=read();
b=read();
n=read();
int t=sqrt(p),flag=0;
ll f=n;//这里最好不要这样写 详见下一份代码
for(int i=1;i<=t;++i)//插入 a^1,a^2,…
{
f=f*b%p;
add(f%10007,f,i);
}
ll g=f=f*qpow(n,p-2,p)%p,k;
for(ll q=t;q-t<p;q+=t,g=g*f%p)//在 a^{√p},a^{2√p},… 中找符合条件的答案
if(k=Find(g%10007,g))
{
printf("%lld\n",q-k);
flag=1;
break;
}
if(!flag)
puts("no solution");
}
return 0;
}
四、exBSGS
如果 $\gcd(a,p)\ne 1$，那么 $a^{k\lfloor\sqrt p\rfloor-r}\equiv b\pmod p$ 就不能把 $a^{-r}$ 移到右边了。
需要想办法把原式转化成底数与模数互质的情况。
根据同余，$a^x\equiv b\pmod p\Leftrightarrow\exists k\in \mathbb Z,a^x+kp=b$。
由于 $\gcd(a,p)\ne 1$，令 $g=\gcd(a,p)$，则 $g|a^x+kp$，当后面一个式子成立时，必然有 $g|b$。
如果此时 $b\nmid a^x+kp$，那么同余方程无解，原方程无解。
否则等式两边同除以 $g$，得到 $\frac{a^x}g+\frac{kp}g=\frac bg$。
由于产生影响的是 $a$，不是 $a^x$，所以此时单独提出来一个 $\frac ag$ 作为“系数”，继续判断是否互质；同理，在 $\frac{kp}g$ 中产生影响的也只有 $p$，所以式子是这样的：$\frac aga^{x-1}+k\frac pg=\frac bg$。
转化回来就是 $\frac aga^{x-1}\equiv\frac bg\pmod{\frac pg}$。又回到了 $a^x\equiv b\pmod p$ 的形式。
一直这样转化，会出现上述无解的情况，或在某一步满足了 $\gcd(a,p)=1$。
这时的式子是
$$
\frac{a^k}{g_1\times g_2\times \cdots\times g_k}a^{x-k}\equiv\frac b{g_1\times g_2\times \cdots\times g_k}\pmod{\frac c{g_1\times g_2\times \cdots\times g_k}}
$$
是 $A\cdot a^{x-k}\equiv B\pmod C$ 的形式。
这样的话直接用 BSGS 的方法就可以解决了。
此外，在这个式子中出现了 $a^{x-k}$，我们无法保证 $x-k\ge 0$。但是如果 $g_i>1$，则 $k\le \lfloor\log_2 C\rfloor$。
这样的话我们可以在程序开始的时候对 $i\in[0,\lfloor\log_2 C\rfloor]$检验 $a^i\equiv b\pmod c$ 是否成立，如果有成立可以直接输出。
这样做还有一个好处。当 $x=0$ 且 $b=1$ 时不受 $\gcd(a,p)$ 是否 $=1$ 的限制。可以直接在检验的时候判掉，不用特判了。注意有的题如果没有定义 $0^0=1$，那么认为 $0^x\equiv 1$ 无解；$0^x\equiv 0$ 的最小非负整数解是 $x=1$。
前面处理到互质的复杂度是 $O(\log n)$，后面是 BSGS，因此总复杂度为 $O(\sqrt n)$。
代码 of POJ 3243
#include<cmath>
#include<cstdio>
#include<cstring>
#define gc() (p1==p2&&(p2=(p1=buf)+fread(buf,1,1<<21,stdin),p1==p2)?EOF:*p1++)
#define ll long long
#define con x=read();z=read();k=read();continue
char buf[2100000],*p1=buf,*p2=buf;
int read()
{
int x=0;
char ch=gc();
while(ch<'0'||ch>'9')
ch=gc();
while(ch>='0'&&ch<='9')
{
x=x*10+(ch&15);
ch=gc();
}
return x;
}
struct edge
{
int x,nxt,v;
}e[100000];
int head[10007],ecnt=-1;
void add(int x,int v)
{
e[++ecnt]=(edge){x,head[x%10007],v};
head[x%10007]=ecnt;
}
ll qpow(ll x,ll y,ll p)
{
ll ans=1;
while(y)
{
if(y&1)
ans=ans*x%p;
x=x*x%p;
y>>=1;
}
return ans;
}
int Find(int x)
{
for(int i=head[x%10007];~i;i=e[i].nxt)
if(e[i].x==x)
return e[i].v;
return 0;
}
int bsgs(int A,int a,int b,int p)
{
for(int i=0;i<10007;++i)
head[i]=-1;
ecnt=-1;
int t=sqrt(p);
ll f=1;
for(int i=1;i<=t;++i)//maybe b doesn't have an inversion
{
f=f*a%p;
add(f*b%p,i);
}
ll g=f*A%p,k;
for(ll q=t;q-t<p;q+=t,g=g*f%p)
if(k=Find(g))
return q-k;
return -1;
}
int gcd(int x,int y)
{
while(y)
{
int t=y;
y=x%y;
x=t;
}
return x;
}
int main()
{
#ifdef wjyyy
freopen("a.in","r",stdin);
#endif
int x,z,k;
x=read();
z=read();
k=read();
while(x|z|k)
{
ll t=1%z;
int flag=0;
for(int i=0;i<=31;++i)
{
if(t==k)
{
printf("%d\n",i);
flag=1;
break;
}
t=t*x%z;
}
if(flag)
{
con;
}
int g=gcd(x,z),cnt=0;
ll A=1;
while(g>1)
{
if(k%g)
{
puts("No Solution");
flag=1;
break;
}
A=A*(x/g)%z;
++cnt;
z/=g;
k/=g;
g=gcd(x,z);
}
if(flag)
{
con;
}
int ans=bsgs(A,x,k,z);
if(ans==-1)
puts("No Solution");
else
printf("%d\n",ans+cnt);
con;
}
return 0;
}
说点什么