[Sdoi2010]古代猪文 （卢卡斯定理，欧拉函数）

哇，这道题真的好好，让我这个菜鸡充分体会到卢卡斯和欧拉函数的强大！

先把题意抽象出来！就是计算这个东西。

p=999911659是素数，p-1=2*3*4679*35617

所以：这样只要求出然后再快速乘法就行了。

那好，怎么做呢？

有模运算的性质得到  然后就是卢卡斯原理。

先把卢卡斯原理放这里：

void init(int mod){                         //对mod取余后，一定小于mod，因此把mod的阶乘存起来就够用
    f[0] = 1;
    for (int i = 1; i <= mod; i++){
        f[i] = f[i - 1] * i % mod;
    }
}



void ex_gcd(LL a, LL b, LL& d, LL& x, LL& y)
{
    if (!b) { d = a; x = 1; y = 0; }
    else{ ex_gcd(b, a%b, d, y, x); y -= x*(a / b); }
}

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



LL Lucas(LL m, LL n, LL p){
    LL  res = 1;
    while (n && m){
        LL n1 = n % p;
        LL m1 = m % p;
        res = res * f[n1] * inv(f[n1 - m1], p) * inv(f[m1], p) % p;
        n /= p;
        m /= p;
    }
    return (res % p + p) % p;
}

则：那么我们其实把它每个存起来Mod[1-4]

然后，就是要找一个值来代替Mod[1-4]。利用中国剩余定理！（哇，太难打了公式了）

什么这样做？因为能同时被2，3，4679，35617那么一定会被99991165同余，那么这个数就是

注意：坑！快速幂一定要加long long，找了3小时的bug

#include<cstdio>
using namespace std;
#define ll long long
const int maxn = 35617;
int N, G, fact[maxn + 5], mod = 999911658;
int prime[5] = { 0, 2, 3, 4679, 35617 }, Mod[5];

void get_fact()

{

    fact[0] = 1;

    for (int i = 1; i <= maxn; i++)

        fact[i] = (ll)fact[i - 1] * i%mod;

}



int ex_t;

void exgcd(int a, int b, int &x, int &y)

{

    if (!b) { x = 1; y = 0; return; }

    exgcd(b, a%b, x, y);

    ex_t = x; x = y; y = ex_t - (a / b)*y;

}



int inv(int a, int p)

{

    int x, y;

    exgcd(a, p, x, y);

    return (x%p + p) % p;

}



int calc(int i, int p)

{

    int ret = 1, x, y, n, m;

    for (x = N, y = i; y; x /= p, y /= p)

    {

        n = x%p; m = y%p;  //卢卡斯定理+预处理阶乘+逆元 

        ret = (ll)ret*fact[n] % p*inv(fact[m], p) % p*inv(n<m ? 0 : fact[n - m], p) % p;

    }

    return ret;

}

ll pow(int x, int n)
{
    int ans = 1;
    for (; n;n>>=1, x=(ll)x*x%mod)
    if (n & 1)ans = (ll)ans*x%mod;
    return ans;
}

int main()
{
    
    scanf("%d%d", &N, &G);
    if (G % (mod + 1) == 0){ printf("0"); return 0; }
    get_fact();            //得到阶乘
    for (int i = 1; i*i <= N; ++i)        //枚举因子
    {
        if (N%i == 0)
        {
            for (int j = 1; j <= 4; ++j)Mod[j] = (Mod[j] + calc(i, prime[j])) % prime[j];
            if (i*i!=N)
            for (int j = 1; j <= 4; ++j)Mod[j] = (Mod[j] + calc(N / i, prime[j])) % prime[j];
        }
    }
    int x, y, b = 0;
    for (int i = 1; i <= 4; i++)  //中国剩余定理 

    {

        exgcd(mod / prime[i], prime[i], x, y);

        b = (b + (ll)Mod[i] % mod*(mod / prime[i]) % mod*x%mod) % mod;

    }
    b = (b%mod + mod) % mod;    mod += 1;
    //printf("%d
", b);
    printf("%lld", pow(G, b));
}