嗯没错,咱又来讲数学了QuQ。
BSGS
原始的BSGS算法是用来解决这样的一类问题:求关于$ x $的方程$$ A^x \equiv B \pmod P $$,其中$ 0 \le x \le P $,且$ P $为质数。
其实算法的思路很简单。
Giant Step
首先记$ m = \lceil \sqrt P \rceil $,则可将$ x $分解为$ x = i \times m + j $,即$ A^x=(A^m)^i \times A^j $,其中$ 0 \le i,j < m $。然后就可以枚举$ i $得到$ j $,$ O(\sqrt P) $级别的枚举。
令$ D=(A^m)^i $,即$ D A^j \equiv B \pmod P $,将$ A^j $视为一个整体解出来,即$ A^j \equiv BD^{-1} \pmod P $。$ D $的逆元可以很方便求出。线性筛逆元那一节我们提到过逆元可积,所以$ A^{-mi}=(A^{-m})^{i} $,所以只要一开始用$ O(\lg P) $计算$ A^{-m} $,接下来直接就有$ D_i^{-1}=D_{i-1}^{-1} \times A^{-m} $。
Baby Step
那得到$ A^j $怎么得到$ j $ 呢?开始预处理的时候只须$ O(\sqrt P) $枚举$ A^j $加入到一个hash表里,然后$ O(1) $查询即可。
总复杂度为$ O(\sqrt P + \lg P) $。
代码(POJ2417,9124K,79ms)
#include <iostream>
#include <cmath>
#include <cstring>
using namespace std;
template<int hashmod = 1 << 20 >
class HashTable
{
private:
int hash[hashmod + 10], stack[1 << 16], stacktop;
int value[hashmod + 10];
int locate(int x)
{
int h = x % hashmod;
while(hash[h] != -1 && hash[h] != x)
h++, h %= hashmod;
return h;
}
public:
HashTable() : stacktop(0){memset(hash, -1, sizeof hash);}
void insert(int x, int v)
{
int h = locate(x);
if(hash[h] == -1)
hash[h] = x, value[h] = v, stack[stacktop++] = h;
}
int get(int x)
{
int h = locate(x);
if(hash[h] != -1)
return value[h];
return -1;
}
void clear(){while(stacktop)hash[stack[--stacktop]] = -1;}
};
HashTable<> hash;
long long power(long long a, long long b, long long c)
{
long long ans = 1;
a %= c;
while(b > 0)
{
if(b & 1)
ans *= a, ans %= c;
a *= a, a %= c;
b >>= 1;
}
return ans;
}
long long BSGS(long long A, long long B, long long C)
{
long long sqrtC = ceil(sqrt(C));
long long base = 1;
hash.clear();
for(int i = 0; i < sqrtC; i++, base *= A, base %= C) //Baby Step
hash.insert(base, i);
long long inv = power(base, C - 2, C);
for(int i = 0, invD = 1, j; i < sqrtC; i++, invD = invD * inv % C) //Giant Step
if((j = hash.get(B * invD % C)) != -1)
return i * sqrtC + j;
return -1;
}
int main()
{
long long A, B, C;
while(cin>>C>>A>>B)
{
long long ans = BSGS(A, B, C);
if(ans == -1)
cout<<"no solution"<<endl;
else
cout<<ans<<endl;
}
return 0;
}
拓展BSGS
拓展后的BSGS,可解决模数不为质数的情况,下文标记为$ C $,即$$ A^x \equiv B \pmod C $$
首先有这样一个性质:$$ ad \equiv bd \pmod {cd} \Leftrightarrow a \equiv b \pmod c $$易证,略。
所以我们可以先执行消因子:
while GCD(A, C) != 1
if B mod GCD(A, C) != 0 then
No Solution
B /= GCD(A, C)
C /= GCD(A, C)
D *= A / GCD(A, C)
cnt ++
注意到其实$ A $是有$ x $个,所以才要执行若干次消元。
若$ A, C $不互质,则不存在$ A $不存在关于$ C $的逆元(证明见此),显然以下无解。
然后问题就变成了$ D \cdot A^{j-cnt} \equiv B \pmod C $。即$ x = i \times m + j + cnt $忽略掉$ cnt $后和上述BSGS就一模一样了。
但还是有一个问题,因为我们默认了$ x \ge cnt $。所以直接$ O(\lg n) $枚举小于cnt的情况,即直接验证是否$ A^i \equiv B \pmod C $,特判即可。
总复杂度仍为$ O(\sqrt C \lg C) $。【好像还是哪里有问题】
代码(POJ3243)
#include <iostream>
#include <cmath>
#include <cstring>
using namespace std;
template<class T = int, int hashmod = 3894229>
class HashTable
{
private:
int hash[hashmod + 10], stack[1 << 16], stacktop;
T value[hashmod + 10];
int locate(int x)
{
int h = x % hashmod;
while(hash[h] != -1 && hash[h] != x)
h++;
return h;
}
public:
HashTable() : stacktop(0){memset(hash, -1, sizeof hash);}
void insert(int x, int v)
{
int h = locate(x);
if(hash[h] == -1)
hash[h] = x, value[h] = v, stack[stacktop++] = h;
}
int get(int x)
{
int h = locate(x);
return hash[h] == x ? value[h] : -1;
}
void clear(){while(stacktop) hash[stack[--stacktop]] = -1;}
};
HashTable<> hash;
void gcd(long long a, long long b, long long& x, long long& y, long long& d)
{
if(b == 0)
d = a, x = 1, y = 0;
else
gcd(b, a % b, y, x, d), y -= x * (a / b);
}
long long gcd(long long a, long long b)
{
return b == 0 ? a : gcd(b, a % b);
}
long long BSGS(long long A, long long B, long long C)
{
for(long long tmp = 1, i = 0; i < 32; i++, tmp = tmp * A % C)
if(tmp == B)
return i;
long long cnt = 0, D = 1;
for(int i; (i = gcd(A, C)) != 1; cnt++)
{
if(B % i != 0)
return -1;
B /= i, C /= i;
D = D * A / i % C;
}
int sqrtC = static_cast<int>(ceil(sqrt(C)));
long long base = 1;
hash.clear();
for(int i = 0; i < sqrtC; i++, base = base * A % C)
hash.insert(base, i);
for(long long i = 0; i < sqrtC; i++)
{
long long x, y, d;
gcd(D, C, x, y, d);
long long j = hash.get((x + C) % C * B % C);
if(j != -1)
return i * sqrtC + j + cnt;
D = D * base % C;
}
return -1;
}
int main()
{
long long A, B, C;
while(cin>>A>>C>>B && (A || B || C))
{
B %= C;
long long ans = BSGS(A, B, C);
if(ans == -1)
cout<<"No Solution"<<endl;
else
cout<<ans<<endl;
}
return 0;
}
POJ3243另一份代码
#include <iostream>
#include <cmath>
#include <cstring>
using namespace std;
template<int hashmod = 1 << 21>
class HashMap
{
private:
int hash[hashmod + 10], stack[1 << 22], stacktop;
int value[hashmod + 10];
int locate(int x)
{
int h = x % hashmod;
while(hash[h] != x && hash[h] != -1)
h++, h %= hashmod;
return h;
}
public:
HashMap():stacktop(0){memset(hash, -1, sizeof hash);}
void insert(int x, int v)
{
int h = locate(x);
if(hash[h] == -1)
hash[h] = x, value[h] = v, stack[stacktop++] = h;
}
int get(int x)
{
int h = locate(x);
if(hash[h] == x)
return value[h];
return -1;
}
void clear(){while(stacktop){hash[stack[--stacktop]] = -1;}}
};
HashMap<> hash;
template<class T>
T gcd(T a, T b, T& x, T& y)
{
if(b == 0)
{
x = 1, y = 0; return a;
}
T ans = gcd(b, a % b, y, x);
y -= x * (a / b);
return ans;
}
template<class T>
T BSGS(T A, T B, T C)
{
T n1, n2;
for(T tmp = 1, i = 0; i < 32; i++, tmp = tmp * A % C)
if(tmp == B)
return i;
T cnt = 0, D = 1;
for(T i; (i = gcd(A, C, n1, n2)) != 1; cnt++)
{
if(B % i != 0)
return -1;
B /= i, C /= i;
D = D * A / i % C;
}
T sqrtC = ceil(sqrt(C));
T base = 1;
hash.clear();
for(T i = 0; i < sqrtC; i++, base = base * A % C)
hash.insert(base, i);
T inv, invD;
gcd(base, C, inv, n1); inv = (inv + C) % C;
gcd(D, C, invD, n1); invD = (invD + C) % C;
for(T i = 0, j; i < sqrtC; i++, invD = invD * inv % C)
if((j = hash.get(invD * B % C)) != -1)
return i * sqrtC + j + cnt;
return -1;
}
int main()
{
long long A, B, C;
while(cin>>A>>C>>B && (A || B || C))
{
B %= C;
long long ans = BSGS(A, B, C);
if(ans != -1)
cout<<ans<<endl;
else
cout<<"No Solution"<<endl;
}
return 0;
}
其它运用
很容易就可以猜到的。其实就是计算对数。求$ \log_a b $的值$ x $。变形即得$ a^x=b $。