HDU 4767( china + 矩阵快速幂)

来源:互联网 发布:淘宝上ios充值的原理 编辑:程序博客网 时间:2024/06/10 21:51

本题目的模数可以拆分啊,95041567 = 31*37*41*43*47;

有bell数有 对任意素数 b[ p+n ] = (b[ n ] + b[ n+1 ]) %p; 这个先用快速幂求出每个小素数模下的值,然后用China搞定。

#include <iostream>#include <cstdio>#include <cstring>#include <algorithm>#include <cmath>typedef long long ll;typedef long long LL;using namespace std;#define rep1(i,x,y) for(int i=x;i<=y;i++)#define repd(i,x,y) for(int i=x;i>=y;i--)#define rep(i,n) for(int i=0;i<(int)n;i++)const int MOD = 95041567;const int N = 100;int mod[]={31,37,41,43,47};int c[5][50][50],A[5];ll ans[5][50];void init(){   rep(k,5){      memset(c[k],0,sizeof(c[k]));      c[k][0][0]=1;      for(int i=0;i<50;i++){          c[k][i][0]=c[k][i][i]=1;          for(int j=1;j<i;j++)             c[k][i][j]=(c[k][i-1][j]+c[k][i-1][j-1])%mod[k];      }   }   rep(k,5){       ans[k][0]=ans[k][1]=1;       rep1(i,2,49){          ans[k][i]=0;          rep1(j,0,i-1)             ans[k][i]=(ans[k][i]+c[k][i-1][j]*ans[k][i-1-j])%mod[k];       }   }}void gcd(ll a,ll b,ll& d,ll& x,ll& y){   if(!b){d=a; x=1; y=0;}   else { gcd(b,a%b,d,y,x); y-=x*(a/b); }}LL china(int n,int* a,int* m){   LL M = 1,d,y,x=0;   rep(i,n) M*=m[i];   for(int i=0;i<n;i++){      LL w=M/m[i];      gcd(m[i],w,d,d,y);      x = (x+y*w*a[i])%M;   }   return (x+M)%M;}int n, modd;struct Matrix{    int mat[N][N];    void show(){      rep(i,n) rep(j,n) {cout<<mat[i][j]<<" "; if(j==n-1) cout<<endl; }    }}a;Matrix Matrix_mul(Matrix a, Matrix b){    Matrix ret;    memset(ret.mat, 0, sizeof(ret.mat));    for(int i = 0; i < n; i++)    for(int j = 0; j < n; j++){        if(a.mat[i][j]){            for(int kk = 0; kk < n; kk++)            ret.mat[i][kk] += (a.mat[i][j]*b.mat[j][kk])%modd;        }    }    return ret;}Matrix Matrix_pow(Matrix tt, int nn){    Matrix ret;    Matrix temp = tt;    memset(ret.mat, 0, sizeof(ret.mat));    for(int i = 0; i< n; i++) ret.mat[i][i] = 1;    while(nn){        if(nn&1) ret = Matrix_mul(ret, temp);        temp = Matrix_mul(temp, temp);        nn >>= 1;    }    return ret;}int nn,te[N];int main(){    init();    int T;    scanf("%d",&T);    while(T--){         scanf("%d",&nn);         if(nn < 50){            rep(i,5) A[i]=ans[i][nn];            printf("%d\n",(int)china(5,A,mod));         }         else {             rep(k,5){                 n=modd=mod[k];                 rep(i,mod[k]-1){                     rep(j,n) a.mat[i][j] = 0;                     a.mat[i][i+1]=1;                 }                 rep(i,mod[k]) a.mat[mod[k]-1][i] = 0;                 a.mat[mod[k]-1][0]=a.mat[mod[k]-1][1]=1;                 for(int i=0,j=1;i<n;i++,j++) te[i]=ans[k][j];                 Matrix tt=Matrix_pow(a,nn-mod[k]);                 A[k]=0;                 rep(i,mod[k]) {A[k]=(A[k]+tt.mat[mod[k]-1][i]*te[i])%mod[k];}             }             printf("%d\n",(int)china(5,A,mod));         }    }    return 0;}


0 0
原创粉丝点击