[SDOI2010]古代猪文 题解

题意简述

给定 \(n,g\)，如果 \(b=\sum\limits_{k|n} \binom{n}{k}\)

求出 \(g^b \bmod 999911659\) 的值

题目分析

为了书写的方便，我们设一个常量 \(M=99911659\)

首先，在 \(g=M\) 的情况下，答案显然为 \(0\)，

否则，根据扩展欧拉定理：

\[\begin{equation} a^b \equiv \begin{cases} a^{b \operatorname{mod} \varphi(m)}, &\gcd(a,m)=1,\\ a^{b}, &\gcd(a,m) \neq 1,b < \varphi(m),\\ a^{(b \operatorname{mod} \varphi(m))+\varphi(m)} &\gcd(a,m) \neq 1,b \leq \varphi(m). \end{cases} \pmod m \end{equation} \]

由于 \(g<10^9\)，所以 \(g<2M\)，再因为 \(M\) 是质数，因此 \(\gcd(g,M)=1\)

再回到扩展欧拉定理，得出 \(g^b\equiv g^{b \operatorname{mod} \varphi(M)}\)

同样因为 \(M\) 为质数，所以 \(\varphi(M)=M-1\)

也就可以将答案转化为 \(g^{b \operatorname{mod} (M-1)}\)，

而要求出这个值，最大的难点是求出 \(b \bmod (M-1)\)

那接下来的问题，就是解决：

\[\sum\limits_{k|n} \binom{n}{k} \bmod 999911658 \]

由于看到了组合数，所以就不得不联想到 Lucas定理

由于直接算显然复杂度为 \(O(n\sqrt{n}+\log M)\)，显然TLE

只能用 Lucas 了，但是 \(M-1\) 不是质数，所以可以用扩展Lucas定理，但由于 \(M-1\) 值太大，又需要优化

但真的有必要这么麻烦吗？

\(M-1\) 是一个常量，所以直接将其质因数分解为 \(999911658=2 \times 3 \times 4679 \times 35617\)

然后将 \(b\) 被 \(2,3,4679,35617\) 的值依次用Lucas算出，然后用中国剩余定理合并就行了

最后再用快速幂算出答案

总时间复杂度为 \(O(\sqrt{n}\log M)\)（注意：\(O(\sqrt{n})\) 是暴力枚举 \(k\) 的复杂度）

参考代码

#include<iostream>
#include<cstdio>
#include<cstring>
using namespace std;
const int MAXN=40007;
const int A=999911658;
const int a[5]={0,2,3,4679,35617};
#define int long long
int factor[MAXN],tot,b[5];

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

int inv(int a,int p)
{
    int x,y;
    exgcd(a,p,x,y);
    return (x%p+p)%p;
}

int qpow(int a,int b,int p)
{
    int ret=1;
    while(b)
    {
        if(b&1) ret=ret*a%p;
        b>>=1;
        a=a*a%p;
    }
    return ret;
}

int fact[MAXN];

int C(int a,int b,int p)
{
    if(b<a) return 0;
    return fact[b]*qpow(fact[a]*fact[b-a]%p,p-2,p)%p;
}

int lucas(int a,int b,int p)
{
    if(!a) return 1;
    else return lucas(a/p,b/p,p)*C(a%p,b%p,p)%p;
}

int n,g;

int cal(int mod)
{
    int ret=0;
    for(int i=1;i<mod;i++) fact[i]=fact[i-1]*i%mod;
    for(int i=1;i<=tot;i++) ret=(ret+lucas(factor[i],n,mod))%mod;
    return ret;
}

int ans;

signed main()
{
    fact[0]=1;
    scanf("%lld%lld",&n,&g);
    if(g==999911659)
    {
        putchar('0');
        return 0;
    }
    for(int i=1;i*i<=n;i++)
    {
        if(n%i==0)
        {
            factor[++tot]=i;
            if(i*i!=n) factor[++tot]=n/i;
        }
    }
    for(int i=1;i<=4;i++) b[i]=cal(a[i]);
    for(int i=1;i<=4;i++)
    {
        int m=A/a[i];
        ans=(ans+m*inv(m,a[i])%A*b[i]%A)%A;
    }
    printf("%lld",qpow(g,ans,A+1));
    return 0;
}