【BZOJ】1951: [Sdoi2010]古代猪文-数论-CRT入门

传送门：bzoj1951

题解

此题考到的知识点比较多，首先得看出这是个$CRT$。
由欧拉定理，当$n\perp p$时，$n^{\phi{p}}\equiv n\ (mod\ p)$.
此题即求$G^{\sum{_{d|n}} \binom{n}{d}(mod\ p-1)}\ (mod \ p)$的值。
$n\leq 10^9$，所以$\binom{n}{d}$可用卢卡斯定理求得：
p为素数时，$\binom{n}{m}\equiv \binom{\lfloor \frac np \rfloor }{\lfloor \frac mp \rfloor}\binom{n \ mod \ p}{m\ mod\ p}\ (mod\ p)$
还得线性推一下逆元。
这里的$p-1（999911658）$不是素数，但它可以拆分成$999911658=2*3*4679*35617$的形式。这样我们可以对四个素因数分别计算，最后通过$CRT$求得最终在$p-1$模意义下的幂次，再快速幂即可。
$CRT$求出的解可以证明，具有唯一性。
$CRT$（中国剩余定理）：求得满足
$x\equiv a_1 (mod\ p_1)$
$x\equiv a_2 (mod\ p_2)$
…
$x\equiv a_n (mod\ p_n)$
的最小正整数$x$。($p_1$~$p_n$为两两互不相同的素数)
设$M=\prod _{i=1} ^n p_i$,$M_i=\frac{M}{p_i}$,$t_i\equiv M_i^{-1} \ (mod\ p_i)$
可以得到$x=\sum _{i=1} ^n M_it_ia_i \ (mod\ M)$
注意：
此题有一种特殊情况，因为欧拉定理只在n,p互质情况下成立，当$n\equiv 0(mod \ p)$时，直接输出0即可。（其实求一遍快速幂也不会错，但乘的过程中，默认$0^0=1$，所以会WA一个点，需要特判一下）

---

代码

#include<cstdio>
#include<algorithm>
#include<cmath>
using namespace std;
typedef long long ll;
const int N=35620,mod=999911659;
int l[4][N],iv[4][N],p[4]={2,3,4679,35617},m[4];
int n,G;

inline int fp(int a,int b,int mod)
{if(!a) return 0;
int ret=1;for(;b;b>>=1,a=1ll*a*a%mod) if(b&1) ret=1ll*ret*a%mod;return ret;}

inline int CRT()
{
    int M=mod-1,ans=0,res;
    for(int i=0;i<4;++i){
        res=1ll*M/p[i]*fp(M/p[i],p[i]-2,p[i])%M*m[i]%M;
        ans=(1ll*ans+res)%M;
    }
    return ans;
}

inline int lucas(int x,int y,int id)
{
    if(!y) return 1;
    if(x<p[id]){
      if(y>x) return 0;
      return 1ll*l[id][x]*iv[id][y]%p[id]*iv[id][x-y]%p[id];
    }else return 1ll*lucas(x/p[id],y/p[id],id)*lucas(x%p[id],y%p[id],id)%p[id];
}

int main(){
    int i,j,k;
    scanf("%d%d",&n,&G);G%=mod;
//  if(G==0){printf("0\n");return 0;}
    for(i=0;i<4;++i) iv[i][0]=l[i][0]=iv[i][1]=1;
    for(i=1;i<4;++i)
     for(j=1;j<p[i];++j){
        l[i][j]=1ll*l[i][j-1]*j%p[i];
        if(j>1) iv[i][j]=1ll*(p[i]-p[i]/j)*iv[i][p[i]%j]%p[i];
     }
    for(i=1;i<4;++i)
     for(j=1;j<p[i];++j)
      iv[i][j]=1ll*iv[i][j]*iv[i][j-1]%p[i];
    for(i=1;i*i<=n;++i){
        if(n%i==0){
            m[0]=(m[0]+((n|i)==n))%2;
            for(j=1;j<4;++j) m[j]=(1ll*m[j]+lucas(n,i,j))%p[j];
            if(i*i!=n) {
                m[0]=(m[0]+((n|(n/i))==n))%2;
                for(j=1;j<4;++j) m[j]=(1ll*m[j]+lucas(n,n/i,j))%p[j];
            }
        }
    }
    printf("%d\n",fp(G,CRT(),mod));
}