bzoj 1951 [Sdoi2010]古代猪文 ——数学综合

题目：https://www.lydsy.com/JudgeOnline/problem.php?id=1951

数学综合题。

费马小定理得指数可以%999911658，又发现这个数可以质因数分解。所以分解做完再用中国剩余定理合并。

为什么不能预处理阶乘的逆元？

为什么正常的中国剩余定理会T？非得两两合并？

　　而且两两合并里的 a0+=m0*x 不太明白。

　　PS：现在明白了。新的a是a=a1+m1*x1=a2+m2*x2，a的通解是a1加上任意倍的m1*x1。

需要特判！那些C( )、lucas( )里的判断也要注意。

#include<iostream>
#include<cstdio>
#include<cstring>
#define ll long long
using namespace std;
const ll m[5]={999911658,2,3,4679,35617};//ll
ll n,g,ans,x,y,a[5],jc[5][35620],jcn[5][35620],M[5];
ll pw(ll x,ll k,ll mod)
{
  ll ret=1;while(k){if(k&1)(ret*=x)%=mod;(x*=x)%=mod;k>>=1;}return ret;
}
void init()
{
  for(int u=1;u<=4;u++)M[u]=m[0]/m[u];
  for(int u=1;u<=4;u++)
    {
      jc[u][0]=1;
      for(int i=1;i<m[u];i++)jc[u][i]=jc[u][i-1]*i%m[u];
//      jcn[u][m[u]-1]=pw(jc[u][m[u]-1],m[u]-2,u);        //为什么不能预处理阶乘的逆元？ 
//      for(int i=m[u]-2;i;i--)jcn[u][i]=(jcn[u][i+1]*(i+1))%m[u];
    }
}
ll C(ll i,ll j,int type)
{
  if(i<j)return 0;    //
  return jc[type][i]*pw(jc[type][j]*jc[type][i-j],m[type]-2,m[type])%m[type];    //
//  if(!j)return 1;
//  return jc[type][i]*jcn[type][j]%m[type]*jcn[type][i-j]%m[type];
}
ll lucas(ll i,ll j,int type)
{
  if(!j||!i)return 1;    //
  if(i<m[type]&&j<m[type])return C(i,j,type);
  return lucas(i/m[type],j/m[type],type)*C(i%m[type],j%m[type],type)%m[type];
}
void exgcd(ll a,ll b)
{
  if(!b){x=1;y=0;return;}
  exgcd(b,a%b);
  ll tp=x;x=y;
  y=tp-a/b*y;
}
int main()
{
  init();
  scanf("%lld%lld",&n,&g);
  if(g==m[0]+1){printf("0");return 0;}    //必须判这个！ 
  for(int i=1;i*i<=n;i++) if(n%i==0)    //
    for(int j=1;j<=4;j++)
    {
        (a[j]+=lucas(n,i,j))%=m[j];
        if(i*i!=n)(a[j]+=lucas(n,n/i,j))%=m[j];
    }
//  ll mod=m[0];            //用中国剩余定理合并：会TLE 
//  for(int i=1;i<=4;i++)
//    {
//      exgcd(M[i],m[i]);
//      (ans+=M[i]*x%mod*a[i]%mod)%=mod;
//    }
//  printf("%lld\n",pw(g,ans,mod+1));
  ll m0=m[1],a0=a[1];        //两个两个地合并 
  for(int i=2;i<=4;i++)
  {
      exgcd(m0,m[i]);
      x=(x*(a[i]-a0)%m[i]+m[i])%m[i];    //%m[i]
      a0+=m0*x;    //
      m0*=m[i];
  }
  printf("%lld\n",pw(g,a0,m[0]+1));
  return 0;
}