[CQOI2018]交错序列
$ solution: $
这一题出得真的很好，将原本一道矩阵快速幂硬生生加入组合数的标签，还那么没有违和感，那么让人看不出来。所以做这道题必须先知道（矩阵快速幂及如何构建递推矩阵）（组合数及二项式定理）。
不知道大家有没有做过洛谷的帕秋莉手环及P哥的桶，这道题中不能有相邻的两个1就是我们在构造这个交错序列时不能连续加入两个1，这个如果直接让我们求方案数（不靠虑一的个数）就是矩阵快速幂的板子了（可以自己推递推方程）。但是这1题偏偏把1的个数搭上了，我们发现1的个数是可以达到 $ 10^8 $ 级别的！所以我们不能直接把1得个数当做一个状态，这里有一个巧妙的避开1的个数的递推方程：
因为我们状态转移时，我们其实只需要知道序列末尾是1还是0，而1的个数我们不妨不管（为什么要折磨自己呢？）我们可以用 $ F[i][j][0/1] $ 为前 $ i $ 位，最后一位为 $ 0/1 $ ，满足条件的序列的1的个数的j次方和（注意是j次方和）。然后我们就可以避开1的个数（因为每一次我们如果要在序列末尾加一个1，就相当于所有长度为i的末尾为一的序列都要加上一个1，然后我们只要求 $ (y+1)^j $ 即可！这个可以用二项式定理搞一下）（这样建出来的矩阵只有90*90可以勉强过）（需要卡常）
后效性：我们避开了一的个数，但是我们需要求零得个数啊！这又怎么办呢？我们把题目要我们求的东西转化一下： $ x^ay^b=(n-y)^ay^b=sum_{i=0}^{a}C(a,i)n^i(-y)^{a-i}y^b $ 这样我们就可以用我们之前设的（满足条件的序列的1的个数的j次方和）的这一维状态求解了！
卡常：因为这一题的模数小于 $ 10^8 $ 不是 $ 10^9 $ 我们矩阵乘法的时候可以多加几次再取模（详情看代码），不然真的很难卡常数！
$ code: $
#include<iostream>
#include<cstdio>
#include<iomanip>
#include<algorithm>
#include<cstring>
#include<cstdlib>
#include<ctime>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#define ull unsigned long long
#define ll long long
#define db double
#define inf 0x7fffffff
#define rg register int
using namespace std;
int n,a,b,m,t,f,tot;
int C[93][93],nf[95];
struct su{
ll s[182][182];
inline void operator *(su x){
su y;
for(rg i=0;i<t;++i)
for(rg j=0;j<t;++j)
y.s[i][j]=0;
for(rg i=0;i<t;++i)
for(rg j=0;j<t;++j){
for(rg k=0;k<t;++k)
y.s[i][j]+=s[i][k]*x.s[k][j];
y.s[i][j]%=m;
}
*this=y;
}
}bas,ans;
inline int qr(){
char ch;
while((ch=getchar())<'0'||ch>'9');
int res=ch^48;
while((ch=getchar())>='0'&&ch<='9')
res=res*10+(ch^48);
return res;
}
int main(){
//freopen(".in","r",stdin);
//freopen(".out","w",stdout);
n=qr(); a=qr(); b=qr(); m=qr();
f=a+b+1; t=f<<1; nf[0]=1;
for(rg i=0;i<f;++i){
C[0][i]=1;
bas.s[i][i]=1;
bas.s[f+i][i]=1;
bas.s[0][f+i]=1;
nf[i+1]=(ull)nf[i]*n%m;
}
for(rg i=1;i<f;++i)
for(rg j=i;j<f;++j)
C[i][j]=(C[i][j-1]+C[i-1][j-1])%m,bas.s[i][f+j]=C[i][j];
ans.s[0][0]=1;
while(n){
if(n&1)ans*bas;
bas*bas; n>>=1;
}
for(rg i=0;i<=a;++i){
rg j=a+b-i;
rg x=(ll)C[i][a]*((a-i)&1?-1:1)*nf[i]%m;
rg y=(ans.s[0][j]+ans.s[0][j+f])%m;
tot=(tot+(ll)x*y)%m;
}printf("%d
",(tot+m)%m);
return 0;
}