exBSGS算法学习笔记

BSGS算法(北上广深算法

用于求解 $a^x\equiv b(modp)$ 的非负整数解 $x$，其中 $(a,p)=1$。

算法过程

令 $x=A\lceil \sqrt{p}\rceil -B$，其中 $A,B\le \lceil \sqrt{p}\rceil$。
那么问题也就转化为找到一组 $A,B$，满足 $a^{A\lceil \sqrt{p}\rceil -B}\equiv b(modp)$。
由于 $(a,p)=1$。
所以 $a^{A\lceil \sqrt{p}\rceil }\equiv b{a}^B(modp)$
暴力枚举 $B$，计算 $b{a}^B$ 的所有值，用哈希表记下该值所对应的最大 $B$。然后暴力枚举 $A$，看 $a^{A\lceil \sqrt{p}\rceil }$ 能否对应到一个 $b{a}^B$。那么一个 $x$ 的值就为 $A\lceil \sqrt{p}\rceil -B$。
大概是一个分块的思想，每 $\sqrt{p}$ 个数为一块进行考虑。
这样做可行是因为 $A\lceil \sqrt{p}\rceil -B$ 覆盖了 $a$ 所有的幂次。
第一个出现的 $A\lceil \sqrt{p}\rceil -B$ 即是最小非负整数解。因为随着 $A$ 增大，答案是递增的。
复杂度是 $O(\sqrt{p})$。

---

---

exBSGS算法

当 $a,p$ 不互质时，怎么求解 $x$？
当 $(a,p)=1$ 时，$a$ 关于 $p$ 的逆元存在，因此可以用 $BSGS$ 求解。
所以，我们要想办法让 $a,p$ 互质。
设 $d_1=\gcd (a,p)$，如果 $d_1\nmid b$，那么方程无解，否则我们将方程同时除掉一个 $d_1$，即：
$$\frac{a}{d_1}{a}^{x-1}\equiv \frac{b}{d_1}(mod\frac{p}{d_1})$$
如果 $a,\frac{p}{d_1}$ 仍然不互质，则进行进行类似操作：
设 $d_2=\gcd (a,\frac{p}{d_1})$，如果 $d_2\nmid b$，那么方程无解，否则我们将方程同时除掉一个 $d_2$，即：
$$\frac{a^2}{d_1{d}_2}{a}^{x-2}\equiv \frac{b}{d_1{d}_2}(mod\frac{p}{d_1{d}_2})$$
假设这样进行 $k$ 次操作后，$a,p$ 终于互质！
记 $D=\prod _{i=1}^k{d}_i$，方程的形式为：
$$\frac{a^k}{D}{a}^{x-k}\equiv \frac{b}{D}(mod\frac{p}{D})$$
这时候，$(a,\frac{p}{D})=1$，于是 $(\frac{a^k}{D},\frac{p}{D})=1$，因此 $\frac{a^k}{D}$ 关于 $\frac{p}{D}$ 的逆元存在。
于是我们将它丢到方程右边，然后就可以 $BSGS$ 求解了！！
这样 $BSGS$ 找到的一个最小正整数解，然后加上 $k$，就是可能的答案了。
不过还可能存在 $x\le k$，使得 $a^x\equiv b(modp)$ 的情况，这种情况我们 $O(k)$ 枚举 $x$ 判断一下就可以了。
总的复杂度还是 $O(\sqrt{p})$

---

---

这是两道模板题

BSGS模板题
exBSGS模板题

代码如下

#include <bits/stdc++.h>
#include<ext/pb_ds/hash_policy.hpp>
#include<ext/pb_ds/assoc_container.hpp>
using namespace __gnu_pbds;
using namespace std;
typedef long long LL;
typedef unsigned long long uLL;
struct custom_hash {
	   static uint64_t splitmix64(uint64_t x) {
	   	x += 0x9e3779b97f4a7c15;
		x = (x ^ (x >> 30)) * 0xbf58476d1ce4e5b9;
		x = (x ^ (x >> 27)) * 0x94d049bb133111eb;
		return x ^ (x >> 31);
	}
	size_t operator()(uint64_t x) const {
		static const uint64_t FIXED_RANDOM = chrono::steady_clock::now().time_since_epoch().count();
		return splitmix64(x + FIXED_RANDOM);
	}
};
LL z = 1;
int read(){
	int x, f = 1;
	char ch;
	while(ch = getchar(), ch < '0' || ch > '9') if(ch == '-') f = -1;
	x = ch - '0';
	while(ch = getchar(), ch >= '0' && ch <= '9') x = x * 10 + ch - 48;
	return x * f;
}
int ksm(int a, int b, int p){
	int s = 1;
	while(b){
		if(b & 1) s = z * s * a % p;
		a = z * a * a % p;
		b >>= 1;
	}
	return s;
}
unordered_map<int, int, custom_hash> f;
void exgcd(int a, int b, int &x, int &y){
	if(!b) x = 1, y = 0;
	else{
		exgcd(b, a % b, y, x);
		y -= a / b * x;
	}
}
int inv(int a, int b){
	int x, y;
	exgcd(a, b, x, y);
	return (x + b) % b;
}
int BSGS(int a, int b, int p){
	int i, j, t, R = sqrt(p - 1) + 1;
	j = b;
	for(i = 0; i <= R; i++){
		f[j] = i;
		j = z * j * a % p;
	}//求出 b * a^i 的值所对应的最大 i 
	j = 1, t = ksm(a, R, p);
	for(i = 0; i <= R; i++){
		if(f.count(j) && i * R > f[j]) return i * R - f[j];//第一次出现即是最小正整数解 
		j = z * j * t % p;
	}
	return -1;//如果无解返回 -1 
}
int exBSGS(int a, int b, int p){
	f.clear();
	int i, ans = -1, g, k = 0, t = 1;
	while(1){
		g = __gcd(a, p);
		if(b % g) return -1;
		if(g == 1) break;
		p /= g; 
		b /= g;//每次都除掉 gcd 
		t = z * t * (a / g) % p;//记录 a^k / D 
		k++;
	}
	b = z * b * inv(t, p) % p;
	ans = BSGS(a, b, p);
	if(ans != -1) ans += k;//求得一个解然后加 k 
	t = 1;
	//可能存在小于等于 k 的解 
	for(i = 0; i <= k; i++){
		if(i % p == b){
			if(ans == -1) ans = i;
			else ans = min(ans, i);
		}
	}
	return ans;
}
int main(){
	int i, j, a, b, p;
	while(scanf("%d", &a) != EOF){
		p = read(); b = read();
		if(!a && !b && !p) break;
		j = exBSGS(a, b, p);
		if(j == -1) printf("No Solution\n");
		else printf("%d\n", j);
	}
	return 0;
}