模板题：POJ2417 调这道题调了我一下午，第二天实在查不出错把#ifndef去掉就A了。。。。。

大致思想

BSGS这个算法主要是用来解决这个问题： Ax≡B(mod C)A^x\equiv B ( mod \ C)已知ABC， 求x

而最普通的BSGS求解的是C为质数（其实就是相当与AC互质）的情况

对于这个问题，如果直接暴力枚举x的值的话，根据Fermat小定理可以知道Aϕp≡1(mod p)A^{\phi p} \equiv 1(mod\ p)

那么我们只要再往下枚举的话必定会产生循环

因此枚举的范围就是00~ϕC\phi C，但是如果C特别大的话显然会超时，于是便有了BSGS算法

算法流程

我们设m=⌈C⌉m = \lceil \sqrt C \rceil（可以证明出只需要取到该值）,

令x=i∗m+jx = i * m + j 那么有Aim+j≡B(mod C)A^{im + j} \equiv B(mod\ C)

因此Aim×Aj≡B(mod C)A^{im} \times A^j\equiv B(mod \ C)

我们枚举i的话，就能得到AimA^{im}的值

令D=AimD=A^{im}，则 D×Aj≡(mod C)D \times A^j \equiv (mod \ C)

显然我们能用扩欧将AjA^j求出，然后可以先预处理出所有AjA^j的值，存在hash中

直接在hash中查找其所对应的j值即可（用map当然可以，只是会T。。。）

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
#include <iostream>
#include <cstdio>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <map>
#define LL long long 
#define ll long long 
using namespace std;
LL A, B, C;
struct hash
{
    static const LL Mod = 100007;
    LL a[Mod + 100], v[Mod + 100];
    LL Find(LL x)
    {
        int y = x % Mod;
        while (a[y] != x && a[y] != -1) (y += 1) %= Mod;
        return y;
    }
    void Insert(LL x, LL y)
    {
        LL z = Find(x);
        if (a[z] == -1)
        {
            a[z] = x;
            v[z] = y;
        }
    }
    LL Get(LL x)
    {
        LL y = Find(x);
        return a[y] == x ? v[y] : -1;
    }
    void Clear()
    {
        memset(a, 0xff, sizeof(a));
    }
} S;
inline LL Pow(LL a, LL b)
{
    LL Ans = 1;
    while (b)
    {
        if (b & 1) (Ans *= a) %= C;
        (a *= a) %= C;
        b >>= 1;
    }
    return Ans;
}
inline LL Ex_gcd(LL a, LL b, LL &x, LL &y)
{
    if (!b)
    {
        x = 1; y = 0;
        return a;
    }
    LL Gcd = Ex_gcd(b, a % b, x, y);
    LL Tmp = x;
    x = y;
    y = (Tmp - (a / b) * y);
    return Gcd;
}

inline int BSGS()
{
    S.Clear();
    LL M = ceil(sqrt((double)C));
    LL Base = 1;
    for (LL i = 0; i < M; ++i)
    {
        S.Insert(Base, i);
        Base = Base * A % C;
    }
    LL D = 1;
    for (LL i = 0; i < M; ++i)
    {
        LL x, y;
        Ex_gcd(D, C, x, y);
        x = ((x * B) % C + C) % C;
        LL k = S.Get(x);
        if (k != -1)
            return i * M + k;
        (D *= Base) %= C;
    }
    return -1;
}
int main()
{
/*#ifndef ONLINE_JDUGE
    freopen("A.in", "r", stdin);
    freopen("A.out", "w", stdout);
#endif
*/
    while (scanf("%lld%lld%lld", &C, &A, &B)!= EOF)
    {
    //  cout<<C<<" "<<A<<" "<<B<<endl;
        LL Ans = BSGS();
        if (Ans == -1) puts("no solution");
        else printf("%lld\n", Ans);
    }
    return 0;
}

扩展

如果A和C不互质怎么办？

设d=gcd(A,C)d = \gcd(A, C)，则A=a∗d,B=b∗d,C=c∗dA=a * d, B = b * d, C = c* d(如果B不是d的倍数且B!=1的话显然无解)

也就是说(a∗d)x≡(b∗d)(mod c∗d)(a * d)^x \equiv (b*d) (mod \ c*d)

因此a∗(a∗d)x−1≡b(mod c)a * (a * d)^{x - 1} \equiv b(mod\ c)

我们把之前的D∗=aD*= a，那么前半部分就只剩下(a∗d)x−1(a * d)^{x - 1} 了

所以我们不断求d=gcd(A,C)d = gcd(A, C),然后不断将B/=d,C/=d,D∗=A/dB /= d, C /= d, D *= A/d，++cnt, 直到d==1为止

最后方程就变成这种形式了：

D∗Ax−cnt≡B(mod C)D * A^{x - cnt} \equiv B(mod\ C)

我们像之前一样令x=i∗m+j+cntx = i * m + j + cnt，和之前一样处理即可

只是这样求出来的x≥cntx \ge cnt，我们需要首先枚举一下0∼cnt−10\sim cnt-1的情况就能解决了