HDU4767 Bell 中国剩余定理贝尔数第二类斯特灵数

首先要了解第二类斯特灵数，还有贝尔数，

接下来是解题的关键：

若p为任意质数，则有Bell(n+p) = Bell(n) + Bell(n+1) (mod p)；

我们来分解题目给的模的数 95041567=31*37*41*43*47，这五个数都是素数，所以可以用上述的式子来解题目，同时又因为这五个数两两互素，所以可以利用中国剩余定理来解 x=a1(mod m1)这个同余方程组，求出 mod（m1*m2...*mi）d的值，这样就可以得到答案了

但是在使用中国剩余定理之前我们要求出 Bell（n）%各个素数的值，这里需要矩阵快速幂

举个例子若素数为5此时可以构造这样一个矩阵

0 1 0 0 0

0 1 1 0 0

0 0 1 1 0

0 0 0 1 1

1 0 0 0 1

再列举一个3的矩阵，你就会发现一定的规律，这样就可以求出 Bell(n)%各个素数的值了

#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
#include  
  
#define ll long long  
#define LL __int64  
#define eps 1e-8  
  
//const ll INF=9999999999999;  
  
#define inf 0xfffffff  
  
using namespace std;  
  
//vector > G;  
//typedef pair P;  
//vector> ::iterator iter;  
//  
//mapmp;  
//map::iterator p;  
//  
//vectorG[30012];  
  
int MOD[5]={31,37,41,43,47};  
  
int a[12][62][62],b[12][62];  
int c[62];  
  
int n;  
  
void getstrling()//获取第二类斯特灵数，可以看维基百科  
{  
    for(int i=0;i<5;i++)  
        a[i][0][0]=1;  
    for(int i=1;i<=50;i++)  
    {  
        for(int j=0;j<5;j++)  
            a[j][i][1]=1;  
        for(int j=1;j<=i;j++)  
        {  
            for(int k=0;k<5;k++)  
                a[k][i][j]=(a[k][i-1][j-1]+j*a[k][i-1][j])%MOD[k];  
        }  
        for(int j=0;j<5;j++)  
        {  
            b[j][i]=0;  
            for(int k=1;k<=i;k++)  
                b[j][i]=(b[j][i]+a[j][i][k])%MOD[j];  
        }  
    }  
}  
  
int Matrixquick(int id,int mod)  
{  
    int x[62],xx[62],y[62][62],yy[62][62];  
    int len=MOD[id]+1;  
    for(int i=1;i<=len+1;i++)  
        for(int j=1;j<=len+1;j++)  
            y[i][j]=0;  
    for(int i=1;i




=n)  
        return x[n];  
    int tempn=n-1;  
    while(tempn)//再写一个函数反而觉得麻烦 直接讲快速幂套在里面了  
    {  
        if(tempn&1)  
        {  
            for(int i=1;i<=len;i++)  
            {  
                xx[i]=0;  
                for(int j=1;j<=len;j++)  
                    xx[i]+=x[j]*y[j][i];  
            }  
            for(int i=1;i<=len;i++)  
                x[i]=xx[i]%mod;  
        }  
        for(int i=1;i<=len;i++)  
        {  
            for(int j=1;j<=len;j++)  
            {  
                yy[i][j]=0;  
                for(int k=1;k<=len;k++)  
                    yy[i][j]+=y[i][k]*y[k][j];  
                yy[i][j]%=mod;  
            }  
        }  
        for(int i=1;i<=len;i++)  
            for(int j=1;j<=len;j++)  
                y[i][j]=yy[i][j];  
        tempn>>=1;  
    }  
    return x[1];  
}  
  
int exgcd(int a,int &x,int b,int &y)  
{  
    if(b==0)  
    {  
        x=1;  
        y=0;  
        return a;  
    }  
    int r=exgcd(b,x,a%b,y);  
    int tmp=x;  
    x=y;  
    y=tmp-a/b*y;  
    return r;  
}  
  
int china()//中国剩余定理  
{  
    int M=1;  
    int ans=0;  
    for(int i=0;i<5;i++)  
        M*=MOD[i];  
    for(int i=0;i<5;i++)  
    {  
        int Mi=M/MOD[i];  
        int x0,y0;  
        int gcd=exgcd(Mi,x0,MOD[i],y0);  
        ans=(ans+x0*Mi*c[i])%M;  
    }  
    return (ans+M)%M;  
}  
  
int main(void)  
{  
    getstrling();  
    int t;  
    scanf("%d",&t);  
    while(t--)  
    {  
        scanf("%d",&n);  
        for(int i=0;i<5;i++)  
            c[i]=Matrixquick(i,MOD[i]);  
        int ans=china();  
        printf("%d\n",ans);  
    }  
}

HDU4767 Bell 中国剩余定理 贝尔数 第二类斯特灵数