luogu_2480: 古代猪文

洛谷：2480古代猪文

题意描述：

给定两个整数(N,G)，求$G^{{sum_{k|n}C_n}k} mod 999911659 $。

数据范围：

(1leq Nleq 10^9,1leq Gleq 10^9)。

思路：

对于这样一个式子，暴力肯定是不可能的，所以我们先来挖掘一些性质。
模数(999911659)是一个质数，我们可以想到对这个式子进行欧拉降幂。
我们可以得到式子：
- (G^{sum_{k|n}C_n^k}equiv G^{sum_{k|n}C_n^kmodvarphi(999911659)} mod(999911659))。
只要当我们求出来指数是多少的时候，那么此时整个结果就变得十分明了了，我们只需要快速幂一下就可以得出结果。
所以接下来我们的重心在于解(sum_{k|n}C_n^kmod varphi(mod))。
对于这样一个式子，我们可以看看(varphi(mod))是多少。
因为(mod=999911659)是一个质数，那么他对应的(varphi(mod)=999911658)。
那么我们这时候的解就变成了解(sum_{k|n}C_n^kmod varphi(999911658))。
因为(n)很大，所以不管怎么算都是很不划算的，我们手里能用的工具只有(lucas)定理，但是(lucas)定理要求模数为质数，这时候我们研究一下(999911658)。
我们发现他虽然不是质数，但是对其质因数分解后，可以发现(999911658=2*3*4679*35617)，四个质数相乘。(出题人设计的真的是太巧妙了。)
于是将上式改写为(sum_{k|n}C_n^kmod varphi(2*3*4679*35617))。
到了这里，似乎就有解头了，因为我们得到了质数，同时对于组合数可以用(lucas)定理来求解。但接下来要怎么考虑呢。因为毕竟上式你不能拆成(原式\%2*原式\%3,...)。
于是我们可以这么考虑，首先设上式为(x)。
- (xequiv num mod(2*3*4679*35617))，其中(num)为上式的最小正整数解。
对于这样一个式子，我们可以拆分成：
- (xequiv num mod(2))
- (xequiv num mod(3))
- (xequiv num mod(4679))
- (xequiv num mod(35617))
可以假设:
- (num\%2=a_1)
- (num\%3=a_2)
- (num\%4679=a_3)
- (num\%35617=a_4)
又因为(num)与(x)同余，所以有：
- (x\%2=a_1)
- (x\%3=a_2)
- (x\%4679=a_3)
- (x\%35617=a_4)
对于(a_i)可以利用(lucas)快速求解，对于(x)可以使用中国剩余定理求解。
假设中国剩余定理求得的结果为(ans)，最终结果就是：
- (G^{ans} mod 999911659)。

注意事项：

可能到了这里就觉得这题能完美解决了，但是其实还有一个小的点需要我们注意，那就是虽然(mod=999911659)是一个质数，但是(g)的范围很大，有可能(g)是(999911659)的倍数，那么不互质就不可以进行常规的欧拉降幂了。
那这里提供两种解决方案。
首先可以一进来就特判一下，看看(g\%999911659==0)否，我们可以特判掉这种情况。
第二种方案可以尝试采用拓展欧拉定理，先来复习一下：
- 当(a,n)不互质时：
  - 当(bgeq varphi(n))，有(a^bequiv a^{bmodvarphi(n)+varphi(n)})。
  - 当(b<varphi(n))，有(a^bequiv a^{bmodvarphi(n)})。
- 那这时候情况很尴尬啊，我们不知道上式的那个组合数求和有多大。
- 但是对于本题，我们可以无脑的采用第一种情况。即(a^bequiv a^{bmodvarphi(n)+varphi(n)})。
- 因为首先当(b)非常大且超过(varphi(n))时，那么式子显然成立。
- 但是当不那么大的时候，对于原式的(mod==999911659)，其实应该采用第二种情况，但是因为(mod)是一个质数，根据欧拉定理(a^{varphi(n)}equiv1mod(n))，其实加上去也就相当于是乘了一个 (1)，对结果没有影响。

代码：

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int mod = 999911659; //是一个质数
ll g, n, prime[10], a[10];

/*
999911658分解质因数结果为 2,3,4679,35617
*/

ll factor[1600], cnt;
void divide(ll x)
{
    for(ll i = 1; i <= x / i; i++)
    {
        if(x % i == 0)
        {
            factor[++cnt] = i;
            if(i != n / i) factor[++cnt] = n / i;
        }
    }
}


ll qmi(ll a, ll b, ll p)
{
    ll res = 1; res %= p;
    while(b)
    {
        if(b & 1) res = res * a % p, res %= p;
        a = a % p * a % p;
        a %= p;
        b >>= 1;
    }
    return res % p;
}

ll C(ll a, ll b, ll p)
{
    if(b > a) return 0;
    if(b > a - b) b = a - b;
    ll x = 1, y = 1;
    for(int i = 0; i < b; i++)
    {
        x = x * (a - i) % p;
        y = y * (i + 1) % p;
    }
    return x * qmi(y, p-2, p);
}

ll lucas(ll a, ll b, ll p)
{
    if(b == 0) return 1;
    return C(a%p, b%p, p) * lucas(a/p, b/p, p) % p;
}

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

ll crt()
{
    ll res = 0;
    for(int i = 1; i <= 4; i++)
        res = (res + a[i] * (mod-1)/prime[i] % (mod-1) * qmi((mod-1)/prime[i], prime[i]-2, prime[i])) % (mod-1);
    return res;
}

int main()
{
    cin >> n >> g;
    /*
    if(g % mod == 0)
    {
        puts("0");
        return 0;
    }
    */
    divide(n);
    memset(a, 0, sizeof a);
    memset(prime, 0, sizeof prime);
    prime[1] = 2; prime[2] = 3;
    prime[3] = 4679; prime[4] = 35617;

    for(int i = 1; i <= 4; i++)
    {
        ll tmp = 0;
        for(int j = 1; j <= cnt; j++)
        {
            ll x = factor[j];
            tmp += lucas(n, x, prime[i]) % prime[i];
            tmp %= prime[i];
        }
        a[i] = tmp;
    }
    ll x = crt();
    //cout << qmi(g, x, mod) << endl;
    cout << qmi(g, x + mod - 1, mod) << endl;
    return 0;
}