需要解决问题之前，首先要做到POJ2417,我的解决问题的方法：http://blog.csdn.net/wyfcyx_forever/article/details/40538515
如今来看这个问题：Ax≡B(mod C)
已知A,B,C<=10^9,给定A,B,C,求x的最小整数解。
注意这里的A,B,C没有不论什么限制！
那么考虑我们的传统的GSBS算法为何不能解决问题：如果枚举的某个i,我们要利用拓展欧几里得求出存不存在某个A^j（0<=j<m）,使得A^(i*m)*A*j%C=B.
那么令A^j=x,我们其实要求的是一个二元方程的整数解：A^(i*m)x+Cy=B.我们知道有解当且仅当gcd(A,C)|B,然而眼下没有不论什么限制，显然是不一定满足的。
我们要对算法进行一些改动。使得其能够进行上述的处理。详细证明去看AekdyCoin犇的题解，我就是简单讲一下我的理解。
一開始的方程等价于A^x*a+C*b=B(a,b是整数)。如今gcd(A，C)!=1,最好还是令t=gcd(A,C)
我们考虑方程(A/t)A^x’%(C/t)=(B/t)的解x’与原来的解x有什么关系。
显然如今的方程等价于(A/t)A^x’*a’+(C/t)*b’=B/t.
两端均乘以t得到：A^(x’+1)*a’+C*b’=B
由系数相等有:x=x’+1,a=a’,b=b’.
那么我们就有一种方法算出x了。先算出x’,再加上1即可了。
考虑怎样求出x’,如今的方程假设A,C依然不互质，就继续迭代将系数除以最大公约数。同一时候前面的常数增大。
然后求出解之后再加上总共除的次数就好了。
Code:
#include <cmath>
#include <cstdio>
#include <cctype>
#include <cstring>
#include <climits>
#include <iostream>
#include <algorithm>
using namespace std;typedef long long LL;
inline LL gcd(LL a, LL b) {
return (!b) ? a : gcd(b, a % b);
}
inline void Exgcd(LL a, LL b, LL &d, LL &x, LL &y) {
if (!b) { d = a, x = 1, y = 0; }
else { Exgcd(b, a % b, d, y, x), y -= x * (a / b); }
}
inline LL Solve(LL a, LL b, LL c) {// ax%c=b S.T. (a,c)=1
LL d, x, y;
Exgcd(a, c, d, x, y);
x = (x + c) % c;
return x * b % c;
}
inline LL Ksm(LL x, LL y, LL p) {
LL res = 1, t = x;
for(; y; y >>= 1) {
if (y & 1) res = res * t % p;
t = t * t % p;
}
return res;
}#define mod 1313131
struct Hashset {
int head[mod], next[35010], f[35010], v[35010], ind;
void reset() {
ind = 0;
memset(head, -1, sizeof head);
}
void Insert(int x, int _v) {
int ins = x % mod;
for(int j = head[ins]; j != -1; j = next[j])
if (f[j] == x) {
v[j] = min(v[j], _v);
return;
}
f[ind] = x, v[ind] = _v;
next[ind] = head[ins], head[ins] = ind++;
}
int operator [] (const int &x) const {
int ins = x % mod;
for(int j = head[ins]; j != -1; j = next[j])
if (f[j] == x)
return v[j];
return -1;
}
}S;LL BSGS(LL C, LL A, LL B, LL p) {// A^x%p=B S.T.(A,p)=1
if (p <= 100) {
LL d = 1;
for(int i = 0; i < p; ++i) {
if (d == B)
return i;
d = d * A % p;
}
return -1;
}
else {
int m = (int)sqrt(p);
S.reset();
LL d = 1, Search;
for(int i = 0; i < m; ++i) {
S.Insert(d, i);
d = d * A % p;
}
for(int i = 0; i * m < p; ++i) {
d = Ksm(A, i * m, p) * C % p;
Search = S[Solve(d, B, p)];
if (Search != -1)
return i * m + Search;
}
return -1;
}
}int main() {
LL x, z, k;
register LL i, j;
while(scanf("%I64d%I64d%I64d", &x, &z, &k) == 3 && (x + z + k)) {
LL d = 1;
bool find = 0;
for(i = 0; i < 100; ++i) {
if (d == k) {
printf("%I64d\n", i);
find = 1;
break;
}
d = d * x % z;
}
if (find)
continue;LL t, C = 1, num = 0;
bool failed = 0;
while((t = gcd(x, z)) != 1) {
if (k % t != 0) {
failed = 1;
break;
}
z /= t;
k /= t;
C = C * x / t % z;
++num;
}if (failed) {
puts("No Solution");
continue;
}LL res = BSGS(C, x, k, z);
if (res == -1)
puts("No Solution");
else
printf("%I64d\n", res + num);
}return 0;
}
版权声明：本文博客原创文章。博客，未经同意，不得转载。
