2023.1.16【数学】BSGS/exBSGS

2023.1.16 [模板]BSGS/exBSGS

全称Boy Step Girl Step

给定一个质数 p，以及一个整数 a，一个整数 b，现在要求你计算一个最小的非负整数 l,

满足$a^x\equiv b(modp)$

算法流程

设t = $\lceil \sqrt{p}\rceil$ ,x = i * t - m (0$\le$i$\le$t , 0$\le$m$<$t - 1)

$a^{i\ast t-b}\equiv b$ (mod p)

$(a^t{)}^i\equiv b\ast a^m$ (mod p)

我们建立一个map,将b * $a^0$ ~ b * $a^{t-1}$ 全部推进去，将值与指数建立映射关系，然后再枚举i，递推算出左式，判断map中是否有值，如果有值则返回答案ans = i * t - map[val]

时间复杂度O($\sqrt{p}$)

inline int BSGS(int a,int b,int p)//a ^ x = b (% p)
{
	map <int,int> q;
	int t = sqrt(p);
	for(int i = 0;i <= t - 1;i++)
		q[(1ll * ksm(a,i,p) * b % p)] = i;
	int c = ksm(a,t,p);
	for(int i = 1;i <= t;i++)
	{
		int now = ksm(c,i,p);
		if(q.find(now) != q.end())
		{
			int x = q[now];
			return i * t - x;
		}
	}
	return -1;
}

exBSGS

那么如果a，p不互质呢？

设d = gcd(a,p);

对于一个同余方程$a^x\equiv b$ (mod p)

我们可以把它转化成二元一次不定方程的形式 $a^x$ + py = b

给a分离系数得到 a$a^{x-1}$ + py = b

想办法消去d，$\frac{a}{d}{a}^{x-1}+\frac{p}{d}y=\frac{b}{d}$

多次迭代，直到d = 1，可以得到以下柿子

$\frac{a^k}{d_1{d}_2...d_k}{a}^{x-k}+\frac{p}{d_1{d}_2...d_k}y=\frac{b}{d_1{d}_2...d_k}$

转化为原式 $\frac{a^k}{d_1{d}_2...d_k}{a}^{x-k}\equiv \frac{b}{d_1{d}_2...d_k}$ $mod(\frac{p}{d_1{d}_2...d_k})$

记录一个变量r = $\frac{a^k}{d_1{d}_2...d_k}$，每次迭代时累乘上$\frac{a}{d_i}$

把r移向同余号右侧，那么我们就得到了一个标准的BSGS，但是同余运算中不能左右同除，于是我们用exgcd算出r关于$\frac{p}{d_1{d}_2...d_k}$的逆元，然后两边同乘以逆元，在跑一边BSGS即可

Code

#include<bits/stdc++.h>
#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/hash_policy.hpp>
typedef long long ll;
using namespace std;
inline int ksm(ll base,int pts,int mod)
{
	ll ans = 1;
	for(;pts > 0;pts >>= 1,base = 1ll * base * base % mod)
		if(pts & 1)
			ans = ans * base % mod;
	return ans;
}
inline int BSGS(int a,ll b,int p)//a ^ x = b (% p)
{
	__gnu_pbds::gp_hash_table <ll,int> q;
	b %= p;
	int t = ceil(sqrt(p));
	for(int i = 0;i <= t;i++)
		q[(ksm(a,i,p) * b % p)] = i;
	ll c = ksm(a,t,p),now = 1;
	for(int i = 1;i <= t;i++)
	{
		now = now * c % p;
		if(q.find(now) != q.end())
		{
			int x = q[now];
			return (i * t - x + p) % p;
		}
	}
	return -1;
}
inline int exgcd(ll &x,int &y,int a,int b)
{
	if(b == 0)
	{
		x = 1;
		y = 0;
		return a;
	}
	int ret = exgcd(x,y,b,a % b);
	int k = y;
	y = x - (a / b) * y;
	x = k;
	return ret;
}
inline int gcd(int a,int b)
{
	if(b == 0) return a;
	return gcd(b,a % b);
}
ll ot = 1;
inline int exBSGS(int a,int b,int p)
{
	if(ot == b) return 0;
	int x,y;
	int d = gcd(a,p);
	if(d == 1)
	{
		ll inv;
		exgcd(inv,y,ot,p);
		inv = ((inv % p) + p) % p;
		return BSGS(a,(1ll * inv * b) % p,p);
	}
	if(b % d != 0) return -1;
	ot = ot * (a / d) % (p / d);
	int ret = exBSGS(a,b / d,p / d);
	if(ret == -1) return -1;
	return ret + 1;
}
signed main()
{
	int a,p,b;
	scanf("%d%d%d",&a,&p,&b);
	while(a || p || b)
	{
		a %= p;b %= p;
		if(b == 1 || p == 1) 
		{
			puts("0");
			scanf("%d%d%d",&a,&p,&b);
			continue; 
		}
		ot = 1;
		int k = exBSGS(a,b,p);
		if(k == -1)
			printf("No Solution\n");
		else
			printf("%lld\n",k);
		scanf("%d%d%d",&a,&p,&b);
	}
	return 0; 
}

语法tip:

C++中有一种比map 和 unordered_map更快的hash映射gp_hash_table，用法和普通map相同

头文件

#include<ext/pb_ds/assoc_container.hpp>
#include<ext/pb_ds/hash_policy.hpp>

命名空间

using namespace __gnu_pbds

定义hash

(__gnu_pbds::)gp_hash_table <int,int> q;