重学BSGS(exBSGS)
重学BSGS(exBSGS)
更好的阅读体验戳此进入
(1 月份写的好像当时忘了发了,现在才发出来。。
目的
BSGS,全称 Baby-Step Giant-Step,即大步小步。一般用于求解离散对数相关问题,即对于同余方程 $ a^x \equiv b \pmod{p} $,给定 $ a, b, p $,求最小的 $ x $,一般要求 $ a \bot p $。其复杂度为 $ O(\sqrt{p}) $。
BSGS
首先我们应该要求 $ x \lt p $,否则可以通过欧拉定理使得 $ x \leftarrow x \bmod \varphi(p) $。
然后我们考虑对 $ x $ 进行拆分,令 $ x = A \lceil \sqrt{p} \rceil - B $,显然 $ 1 \le A, B \le \lceil \sqrt{p} \rceil $。
然后带入原式,即:$ a^{A \lceil \sqrt{p} \rceil - B} \equiv b \pmod{p} $。
将 $ B $ 移项得到:$ a^{A \lceil \sqrt{p} \rceil} \equiv b \times a^B \pmod{p} $
于是我们可以考虑首先枚举所有 $ B $,然后通过 unordered_map
存储所有 $ b \times a^B \bmod p \rightarrow B $ 的映射,然后再枚举 $ A $,对于每个 $ a^{A \lceil \sqrt{p} \rceil} $ 都在映射中尝试寻找对应的 $ B $,若存在,此时则有 $ A, B $ 的取值,直接带入原式即可得到 $ x $ 的取值。
同时此处有个小 trick,我们正序枚举 $ A $ 的时候,因为 $ B \le \lceil \sqrt{p} \rceil $,所以枚举过程中遇到的第一个合法 $ x $ 便一定是最小的合法 $ x $,可以直接输出。
例题 #1
题面
LG-P3846 [TJOI2007] 可爱的质数/【模板】BSGS。
BSGS 模板。
Solution
BSGS 模板。
对于快速幂部分可以使用光速幂以根号预处理并 $ O(1) $ 求解,但本题复杂度可以多一只 $ \log $,故没必要。
Code
#define _USE_MATH_DEFINES
#include <bits/stdc++.h>
#define PI M_PI
#define E M_E
#define npt nullptr
#define SON i->to
#define OPNEW void* operator new(size_t)
#define ROPNEW(arr) void* Edge::operator new(size_t){static Edge* P = arr; return P++;}
using namespace std;
mt19937 rnd(random_device{}());
int rndd(int l, int r){return rnd() % (r - l + 1) + l;}
bool rnddd(int x){return rndd(1, 100) <= x;}
typedef unsigned int uint;
typedef unsigned long long unll;
typedef long long ll;
typedef long double ld;
template < typename T = int >
inline T read(void);
int P, B, N;
unordered_map < int, int > mp;
ll qpow(ll a, ll b){
ll ret(1), mul(a);
while(b){
if(b & 1)ret = ret * mul % P;
b >>= 1;
mul = mul * mul % P;
}return ret;
}
int main(){
P = read(), B = read(), N = read();
int lim = (int)ceil(sqrt(P));
for(int i = 1; i <= lim; ++i)mp.insert({(ll)N * qpow(B, i) % P, i});
for(int i = 1; i <= lim; ++i){
ll val = qpow(B, (ll)i * lim);
if(mp.find(val) != mp.end())printf("%lld\n", (ll)i * lim - mp[val]), exit(0);
}puts("no solution\n");
fprintf(stderr, "Time: %.6lf\n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}
template < typename T >
inline T read(void){
T ret(0);
int flag(1);
char c = getchar();
while(c != '-' && !isdigit(c))c = getchar();
if(c == '-')flag = -1, c = getchar();
while(isdigit(c)){
ret *= 10;
ret += int(c - '0');
c = getchar();
}
ret *= flag;
return ret;
}
例题 #2
题面
给定 $ m, k $,求满足 $ 111 \cdots 111 \equiv k \pmod{m} $(共 $ n $ 个 $ 1 $)的最小 $ n $。保证 $ m $ 为质数。
Solution
尝试进行转化,对 $ \texttt{LHS} \leftarrow \texttt{LHS} \times 9 + 1 $ 即可转换为 $ 10 $ 的幂,即:
然后套一下 BSGS 模板即可。
依然可以不使用光速幂,最终复杂度 $ O(\sqrt{m} \log m) $。
注意过程中部分需要使用 __int128_t
。
Code
#define _USE_MATH_DEFINES
#include <bits/stdc++.h>
#define PI M_PI
#define E M_E
#define npt nullptr
#define SON i->to
#define OPNEW void* operator new(size_t)
#define ROPNEW(arr) void* Edge::operator new(size_t){static Edge* P = arr; return P++;}
using namespace std;
mt19937 rnd(random_device{}());
int rndd(int l, int r){return rnd() % (r - l + 1) + l;}
bool rnddd(int x){return rndd(1, 100) <= x;}
typedef unsigned int uint;
typedef unsigned long long unll;
typedef long long ll;
typedef long double ld;
template < typename T = int >
inline T read(void);
ll K, M;
unordered_map < ll, ll > mp;
ll qpow(ll a, ll b){
ll ret(1), mul(a);
while(b){
if(b & 1)ret = (__int128_t)ret * mul % M;
b >>= 1;
mul = (__int128_t)mul * mul % M;
}return ret;
}
int main(){
K = read < ll >(), M = read < ll >();
ll lim = (ll)ceil(sqrt(M));
for(int i = 1; i <= lim; ++i)mp.insert({(__int128_t)(9 * K + 1) % M * qpow(10, i) % M, i});
for(int i = 1; i <= lim; ++i){
ll val = qpow(10, (ll)i * lim);
if(mp.find(val) != mp.end())printf("%lld\n", (ll)i * lim - mp[val]), exit(0);
}
fprintf(stderr, "Time: %.6lf\n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}
template < typename T >
inline T read(void){
T ret(0);
int flag(1);
char c = getchar();
while(c != '-' && !isdigit(c))c = getchar();
if(c == '-')flag = -1, c = getchar();
while(isdigit(c)){
ret *= 10;
ret += int(c - '0');
c = getchar();
}
ret *= flag;
return ret;
}
例题 #3
题面
给定 $ y, z, p $,求一些东西。
Solution
快速幂,乘法逆元,BSGS。
其中有个地方需要注意,即若 $ y \bot p $ 不成立,那么若 $ z \not\equiv 0 $ 则无解,若 $ z \equiv 0 $ 则解为 $ 1 $。
Code
#define _USE_MATH_DEFINES
#include <bits/stdc++.h>
#define PI M_PI
#define E M_E
#define npt nullptr
#define SON i->to
#define OPNEW void* operator new(size_t)
#define ROPNEW(arr) void* Edge::operator new(size_t){static Edge* P = arr; return P++;}
using namespace std;
mt19937 rnd(random_device{}());
int rndd(int l, int r){return rnd() % (r - l + 1) + l;}
bool rnddd(int x){return rndd(1, 100) <= x;}
typedef unsigned int uint;
typedef unsigned long long unll;
typedef long long ll;
typedef long double ld;
template < typename T = int >
inline T read(void);
ll Y, Z, P;
unordered_map < ll, ll > mp;
ll qpow(ll a, ll b){
ll ret(1), mul(a);
while(b){
if(b & 1)ret = ret * mul % P;
b >>= 1;
mul = mul * mul % P;
}return ret;
}
int main(){
int T = read(), K = read();
while(T--){
Y = read(), Z = read(), P = read();
mp.clear();
switch(K){
case 1:{
printf("%lld\n", qpow(Y, Z));
break;
}
case 2:{
if(Y % P == 0 && Z % P){printf("Orz, I cannot find x!\n"); break;}
printf("%lld\n", Z * qpow(Y, P - 2) % P);
break;
}
case 3:{
if(Y % P == 0 && Z % P){printf("Orz, I cannot find x!\n"); break;}
if(Y % P == 0 && Z % P == 0){printf("1\n"); break;}
ll lim = (ll)ceil(sqrt(P));
for(int i = 1; i <= lim; ++i)mp[Z * qpow(Y, i) % P] = i;
for(int i = 1; i <= lim; ++i){
ll val = qpow(Y, (ll)i * lim);
if(mp.find(val) != mp.end()){printf("%lld\n", (ll)i * lim - mp[val]); break;}
if(i == lim)printf("Orz, I cannot find x!\n");
}
break;
}
default: break;
}
}
fprintf(stderr, "Time: %.6lf\n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}
template < typename T >
inline T read(void){
T ret(0);
int flag(1);
char c = getchar();
while(c != '-' && !isdigit(c))c = getchar();
if(c == '-')flag = -1, c = getchar();
while(isdigit(c)){
ret *= 10;
ret += int(c - '0');
c = getchar();
}
ret *= flag;
return ret;
}
例题 #4
题面
给定 $ p, a, b, x_1, t $,令:
求 $ x_n = t $ 的最小 $ n $。
Solution
将原式通过等比数列求和公式硬推一下即可最终处理成 $ a^n \equiv b \pmod{p} $ 的形式,然后用 BSGS 即可。
但是需要注意这东西细节巨多,如等比数列求和公比是否为 $ 1 $,部份数是否为 $ 0 $ 等,很多都需要分类讨论。
Code
#define _USE_MATH_DEFINES
#include <bits/stdc++.h>
#define PI M_PI
#define E M_E
#define npt nullptr
#define SON i->to
#define OPNEW void* operator new(size_t)
#define ROPNEW(arr) void* Edge::operator new(size_t){static Edge* P = arr; return P++;}
using namespace std;
mt19937 rnd(random_device{}());
int rndd(int l, int r){return rnd() % (r - l + 1) + l;}
bool rnddd(int x){return rndd(1, 100) <= x;}
typedef unsigned int uint;
typedef unsigned long long unll;
typedef long long ll;
typedef long double ld;
template < typename T = int >
inline T read(void);
ll P, A, B, X1, T;
unordered_map < ll, ll > mp;
ll qpow(ll a, ll b){
ll ret(1), mul(a);
while(b){
if(b & 1)ret = ret * mul % P;
b >>= 1;
mul = mul * mul % P;
}return ret;
}
int main(){
int t = read();
while(t--){
mp.clear();
P = read(), A = read(), B = read(), X1 = read(), T = read();
if(X1 == T){printf("1\n"); continue;}
if(A == 0){printf("%d\n", B == T ? 2 : -1); continue;}
if(A == 1){
if(B == 0){printf("-1\n"); continue;}
printf("%lld\n", ((T - X1) % P + P) % P * qpow(B, P - 2) % P + 1);
continue;
}
ll b = ((T - B * qpow((1 - A + P) % P, P - 2) % P) % P + P) % P * qpow(((X1 - B * qpow((1 - A + P) % P, P - 2) % P) % P + P) % P, P - 2) % P;
ll lim = (ll)ceil(sqrt(P));
for(int i = 1; i <= lim; ++i)mp.insert({b * qpow(A, i) % P, i});
for(int i = 1; i <= lim; ++i){
ll val = qpow(A, (ll)i * lim);
if(mp.find(val) != mp.end()){printf("%lld\n", (ll)i * lim - mp[val] + 1); break;}
if(i == lim)puts("-1");
}
}
fprintf(stderr, "Time: %.6lf\n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}
template < typename T >
inline T read(void){
T ret(0);
int flag(1);
char c = getchar();
while(c != '-' && !isdigit(c))c = getchar();
if(c == '-')flag = -1, c = getchar();
while(isdigit(c)){
ret *= 10;
ret += int(c - '0');
c = getchar();
}
ret *= flag;
return ret;
}
exBSGS
与普通 BSGS 的区别就是不保证 $ p $ 为质数,或者说不保证 $ a \bot p $。
大致思路就是先考虑将同余方程转化为普通方程,即:
我们令 $ d = \gcd(a, p) $。发现必须满足 $ d \mid b $,否则无解,则可转化为:
换一种写法:
写成同余方程就是:
此时若满足 $ a \bot \dfrac{p}{d} $ 那么则可以直接进行 BSGS,否则继续递归下去再拆一次,直到两者互质为止,然后再将答案转换回对应的 $ x $ 即可。
然后还有一个问题就是,对于递归层数为 $ m $ 时,我们只能得到 $ \ge m $ 的最小解,所以若最小解在 $ [0, m) $ 的话是不会被统计到的,所以我们还需要枚举验证一下这一部分。
同时注意过程中还需要特判若 $ b = 1 $ 那么直接返回 $ 0 $ 即可。
然后具体的写法可以考虑写成非递归形式,会更简洁,常数也会更小。
例题 #5
题面
exBSGS 模板。
Solution
exBSGS 模板。
Code
#define _USE_MATH_DEFINES
#include <bits/stdc++.h>
#define PI M_PI
#define E M_E
#define npt nullptr
#define SON i->to
#define OPNEW void* operator new(size_t)
#define ROPNEW(arr) void* Edge::operator new(size_t){static Edge* P = arr; return P++;}
using namespace std;
mt19937 rnd(random_device{}());
int rndd(int l, int r){return rnd() % (r - l + 1) + l;}
bool rnddd(int x){return rndd(1, 100) <= x;}
typedef unsigned int uint;
typedef unsigned long long unll;
typedef long long ll;
typedef long double ld;
template < typename T = int >
inline T read(void);
map < ll, ll > mp;
ll A, P, B;
ll BSGS(ll a, ll b, ll MOD, ll k = 1){
mp.clear();
ll lim = (ll)ceil(sqrt(MOD));
ll cur(1);
for(int i = 1; i <= lim; ++i)(cur *= a) %= MOD, mp[b * cur % MOD] = i;
ll powA = cur * k % MOD;
for(int i = 1; i <= lim; ++i){
if(mp.find(powA) != mp.end())return i * lim - mp[powA];
(powA *= cur) %= MOD;
}return -0x3f3f3f3f;
}
ll exBSGS(ll a, ll b, ll MOD){
ll originA = a % MOD, originB = b % MOD, originMOD = MOD;
a %= MOD, b %= MOD;
if(b == 1 || MOD == 1)return 0;
ll cnt(0), k(1);
ll cur(1);
while(true){
if(cur == originB)return cnt;
(cur *= originA) %= originMOD;
ll d = __gcd(a, MOD);
if(b % d)return -0x3f3f3f3f;
if(d == 1)return BSGS(a, b, MOD, k) + cnt;
k = k * a / d % MOD, b /= d, MOD /= d;
++cnt;
}
}
int main(){
while(true){
A = read(), P = read(), B = read();
if(!A && !P && !B)exit(0);
ll ret = exBSGS(A, B, P);
if(ret < 0)printf("No Solution\n");
else printf("%lld\n", ret);
}
fprintf(stderr, "Time: %.6lf\n", (double)clock() / CLOCKS_PER_SEC);
return 0;
}
template < typename T >
inline T read(void){
T ret(0);
int flag(1);
char c = getchar();
while(c != '-' && !isdigit(c))c = getchar();
if(c == '-')flag = -1, c = getchar();
while(isdigit(c)){
ret *= 10;
ret += int(c - '0');
c = getchar();
}
ret *= flag;
return ret;
}
UPD
update-2023_01_06 初稿