AcWing 213. 古代猪文 数学知识

Posted 2022-07-01 l999q

传送门

题目描述：

给定整数n,q，计算 $q^\sum_d|n C_n^d$ mod 999911659。

输入格式

输入包括一行，包含两个整数n，q，用一个空格隔开。

输出格式

输出包括一行，包含一个整数表示最终结果。

数据范围

1≤n,q≤10⁹ 

输入样例：

4 2

输出样例：

2048

 提示:对于n的每一个正因数d，都有一个$C_n^d$的值，将它们全部加起来得到的和就是$\sum_d|n C_n^d$ 。

 

题解：经典题，用到的数学知识比较多。

　　当q为999911659时，答案为0。否则，因为999911659是质数，由欧拉定理推论得：

 

　　$q^\sum_d|n C_n^d \equiv q^\sum_d|n C_n^d mod 999911658 mod 999911659$

 

　　因此，本题的关键是计算

 

　　$\sum_d|n C_n^d mod 999911658$。

 

　　n很大，分解质因子可以发现999911658 = 2 * 3 * 4679 * 35617，我们可以用Lucas定理求组合数 $C_n^d$ ，分别计算出

 

　　$\sum_d|n C_n^d$ 对2,3,4679,35617四个质数取模的结果，记为d[1],d[2],d[3],d[4]。

 

　　最后用中国剩余定理求解线性同余方程。

 

 $\left\\beginmatrix x & mod & 2 & = & a_1\\  x & mod & 3 & = & a_2 \\  x & mod & 4679 & = & a_3 \\  x & mod & 35617 & = & a_4 \endmatrix\right.$

 

代码： 

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int N = 4e4 + 10;
const ll mod = 999911659;
ll fac[N] = 1,n,q,b[5];
int a[5] = 0,2,3,4679,35617;
ll qp(ll a,ll b) 
    ll ans = 1;
    for (;b;b>>=1,a=a*a%mod) 
        if (b&1) ans = ans * a %mod;
    return ans;

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

ll getinv(ll a,ll b) 
    ll x,y;
    exgcd(a,b,x,y);
    return (x%b+b)%b;

ll lucas(ll x,ll p) 
    ll ans = 1;
    for (ll y = n;x;x/=p,y/=p) 
        ll a = x%p, b = y%p;
        ans = ans*fac[b]%p*getinv(fac[a],p)%p*getinv(b<a?0:fac[b-a],p)%p;
    
    return ans;

int main()
    scanf("%lld%lld",&n,&q);
    q%=mod;
    if(!q) 
        printf("0\n");
        return 0;
    
    ll mod1 = mod -1,ans = 0,x,y;
    for (ll i = 1; i <= 35617; i++) fac[i] = fac[i-1]*i%mod1;
    for (ll i = 1; i*i<=n; i++)
        if (n%i==0) 
            for (int j = 1; j <= 4; j++)
                b[j] = (b[j]+lucas(i,a[j]))%a[j];
            if (i*i!=n) for (int j = 1; j <= 4; j++)
                b[j] = (b[j]+lucas(n/i,a[j]))%a[j];
        
    for (int i = 1; i <= 4; i++)
        exgcd(mod1/a[i],a[i],x,y);
        ans = (ans+x*(mod1/a[i])%mod*b[i])%mod1;
    
    ans = (ans+mod1)%mod1;
    ans = qp(q,ans);
    printf("%lld\n", ans);
    return 0;

View Code