[BZOJ5298][CQOI2018]交错序列(DP+矩阵乘法)
https://blog.csdn.net/dream_maker_yk/article/details/80377490
斯特林数有时并没有用。
1 #include<cstdio> 2 #include<cstring> 3 #include<algorithm> 4 #define rep(i,l,r) for (int i=(l); i<=(r); i++) 5 typedef long long ll; 6 using namespace std; 7 8 const int N=210; 9 int n,a,b,mod,m,k,ans,fac[N],inv[N]; 10 11 struct Mat{ 12 int a[N][N]; 13 Mat(){ memset(a,0,sizeof(a)); } 14 }; 15 16 Mat operator *(const Mat &a,const Mat &b){ 17 Mat c; 18 rep(i,0,m-1) rep(j,0,m-1) if (a.a[i][j]) 19 rep(k,0,m-1) if (b.a[j][k]) c.a[i][k]=(c.a[i][k]+1ll*a.a[i][j]*b.a[j][k])%mod; 20 return c; 21 } 22 23 Mat ksm(const Mat &a,int b){ 24 Mat c=a,res; 25 rep(i,0,m-1) res.a[i][i]=1; 26 for (; b; c=c*c,b>>=1) 27 if (b & 1) res=res*c; 28 return res; 29 } 30 31 int C(int n,int m){ return n<m ? 0 : 1ll*fac[n]*inv[m]%mod*inv[n-m]%mod; } 32 33 int ksm(int a,int b){ 34 int res=1; 35 for (; b; a=1ll*a*a%mod,b>>=1) 36 if (b & 1) res=1ll*res*a%mod; 37 return res; 38 } 39 40 int main(){ 41 scanf("%d%d%d%d",&n,&a,&b,&mod); 42 k=a+b+1; m=2*k; 43 fac[0]=1; rep(i,1,m) fac[i]=1ll*fac[i-1]*i%mod; 44 inv[m]=ksm(fac[m],mod-2); 45 for (int i=m-1; ~i; i--) inv[i]=1ll*inv[i+1]*(i+1)%mod; 46 Mat s; rep(i,0,k-1) s.a[i][i+k]=1; 47 rep(i,0,k-1) rep(j,0,i) s.a[j][i]=s.a[j+k][i]=C(i,j); 48 s=ksm(s,n); int x=1,ans=0; 49 rep(i,0,b) ans=(ans+1ll*C(b,i)*x%mod*(s.a[0][a+b-i]+s.a[0][a+b-i+k])%mod*(((b-i)&1)?-1:1))%mod,x=1ll*x*n%mod; 50 printf("%d\n",(ans+mod)%mod); 51 return 0; 52 }
