SDOI2010：古代猪文（Lucas定理 + CRT合并）

 

题目描述

猪王国的文明源远流长，博大精深。

iPig在大肥猪学校图书馆中查阅资料，得知远古时期猪文文字总个数为N。当然，一种语言如果字数很多，字典也相应会很大。当时的猪王国国王考虑到如果修一本字典，规模有可能远远超过康熙字典，花费的猪力、物力将难以估量。故考虑再三没有进行这一项劳猪伤财之举。当然，猪王国的文字后来随着历史变迁逐渐进行了简化，去掉了一些不常用的字。

iPig打算研究古时某个朝代的猪文文字。根据相关文献记载，那个朝代流传的猪文文字恰好为远古时期的k分之一，其中k是N的一个正约数（可以是1和N）。不过具体是哪k分之一，以及k是多少，由于历史过于久远，已经无从考证了。

iPig觉得只要符合文献，每一种能整除N的k都是有可能的。他打算考虑到所有可能的k。显然当k等于某个定值时，该朝的猪文文字个数为N / k。然而从N个文字中保留下N / k个的情况也是相当多的。iPig预计，如果所有可能的k的所有情况数加起来为P的话，那么他研究古代文字的代价将会是G的P次方。

现在他想知道猪王国研究古代文字的代价是多少。由于iPig觉得这个数字可能是天文数字，所以你只需要告诉他答案除以999911659的余数就可以了。

输入输出格式

输入格式：

 

 

 

输入文件ancient.in有且仅有一行：两个数N、G，用一个空格分开。

 

输出格式：

 

 

 

输出文件ancient.out有且仅有一行：一个数，表示答案除以999911659的余数。

 

输入输出样例

输入样例#1： 复制

4 2

输出样例#1： 复制

2048

说明

数据规模

10%的数据中，1 <= N <= 50；

20%的数据中，1 <= N <= 1000；

40%的数据中，1 <= N <= 100000；

100%的数据中，1 <= G <= 1000000000，1 <= N <= 1000000000。

 

题意：计算

 

思路：根据费马小定理，计算指数部分时取模999911658即可（特判G=999911659，此时不适用费马小定理），由于数据范围很大，留意到999911658=2*3*4679*35617，因此对每个质因数用Lucas定理计算，最后用中国剩余定理合并结果，又因为模数均为质因子，CRT求逆元时可以用费马小定理计算。

# include <bits/stdc++.h>
using namespace std;
typedef long long LL;
const int MOD = 999911659, mod = 999911658;
int prime[] = {2, 3, 4679, 35617};
int fac[40000];//比模数要大
int num[4];
int qmod(LL x, int y, int p)
{
    int ans = 1;
    for(;y;y>>=1)
    {
        if(y&1) ans = ans*x%p;
        x = x*x%p;
    }
    return ans;
}
int C(int n, int m, int p)
{
    if (n < m) return 0;
    return fac[n]*qmod(fac[n-m]*fac[m], p-2, p)%p;
}

int lucas(int n, int m, int p)
{
    if (m == 0) return 1;
    return C(n%p, m%p, p)*lucas(n/p, m/p, p)%p;
}
int crt()
{
    int ans = 0;
    for(int i=0; i<4; i++)
    {
        int Mi = mod/prime[i];
        ans = (ans + 1LL* Mi * qmod(Mi, prime[i]-2, prime[i])%mod * num[i] %mod) % mod;
    }
    return ans;
}
int main()
{
    int n, g;
    scanf("%d%d",&n,&g);
    if(g == MOD) return 0*puts("0");
    for(int i=0; i<4; ++i)
    {
        fac[0] = 1;
        for (int j = 1; j <= prime[i]; ++j)
            fac[j] = fac[j-1]*j%prime[i];
        for(int j=1;j<=sqrt(n); ++j)
        {
            if(n%j == 0)
            {
                num[i] = (num[i]+lucas(n, j, prime[i]))%prime[i];
                if(j*j!=n) num[i] = (num[i]+lucas(n, n/j, prime[i]))%prime[i];
            }
        }
    }
    printf("%d\n",qmod(g, crt(), MOD));
    return 0;
}