POJ 3243 Clever Y(扩展BSGS,gcd(a,p)!=1)

题目链接：
POJ 3243 Clever Y
题意：
跟POJ 2417 Discrete Logging类似，只不过$gcd(a,p)!=1$.
分析；
初始化$cnt = 0$(消因子轮数),$d = 1$(消掉的$gcd$乘积).
令$tmp = gcd(a, p)$ 当$tmp != 1$时,修改变量值:$b /= tmp$(先判断$b$是否是$tmp$的倍数),$p /= tmp$, $d = a/tmp * d % p;$
通过若干轮消掉$a, p$的因子使得最终$gcd(a, p) = 1$.这时再调用普通的$BSGS$得到解为$res$，则最终答案是$res + cnt$。
但是这样求得的解是$\geq cnt$的，需要先判断下是否有 <cnt <script type="math/tex" id="MathJax-Element-64"> $cnt$的最大值。因为每次消去的最小因子是 $2$，可以得到$cnt$的最大值是 $log2(p)$，先跑一遍 $50$次遍历的循环是绰绰有余的。

#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
#include <algorithm>
#include <climits>
#include <cmath>
#include <ctime>
#include <cassert>
#define IOS ios_base::sync_with_stdio(0); cin.tie(0);
using namespace std;
typedef long long ll;

const ll MOD = 100007;
ll hs[MOD + 100], id[MOD + 100];

ll find(ll x)
{
    ll t = x % MOD;
    while(hs[t] != x && hs[t] != -1) t = (t + 1) % MOD;
    return t;
}

void insert(ll x, ll ii)
{
    ll pos = find(x);
    if(hs[pos] == -1){
        hs[pos] = x;
        id[pos] = ii;
    }
}

ll get(ll x)
{
    ll pos = find(x);
    return hs[pos] == x ? id[pos] : -1;
}

ll ex_gcd(ll a, ll b, ll& x, ll& y)
{
    if(b == 0) {
        x = 1, y = 0;
        return a;
    }
    ll d = ex_gcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}

ll inv(ll a, ll p)
{
    ll d, x, y;
    d = ex_gcd(a, p, x, y);
    return d == 1 ? (x % p + p) % p : -1;
}

ll BSGS(ll a, ll b, ll p)
{
    memset(hs, -1, sizeof(hs));
    memset(id, -1, sizeof(id));
    ll m = (ll)ceil(sqrt(p + 0.5));
    ll tmp = 1;
    for(ll i = 0; i < m; i++) {
        insert(tmp, i);
        tmp = tmp * a % p;
    }
    ll base = inv(tmp, p);
    ll res = b, z;
    for(ll i = 0; i < m; i++) {
        if((z = get(res)) != -1) return i * m + z;
        res = res * base % p;
    }
    return -1;
}

ll gcd(ll x, ll y)
{
    return y == 0 ? x : gcd(y, x % y);
}

ll solve(ll a, ll b, ll p)
{
    ll tmp = 1;
    for(int i = 0; i <= 50; ++ i) {
        if(tmp == b) return i;
        tmp = tmp * a % p;
    }
    ll cnt = 0, d = 1 % p;
    while((tmp = gcd(a, p)) != 1) {
        if(b % tmp) return -1;
        b /= tmp;
        p /= tmp;
        d = a / tmp * d % p;
        cnt++;
    }
    b = b * inv(d, p) % p;
    ll ans = BSGS(a, b, p);
    if(ans == -1) return -1;
    else return ans + cnt;
}

int main()
{
    ll a, b, p;
    while(~scanf("%lld%lld%lld", &a, &p, &b) && (a || p || b)){
        ll ans = solve(a, b % p, p);
        if(ans == -1) printf("No Solution\n");
        else printf("%lld\n", ans);
    }
    return 0;
}