[CQOI2018]交错序列

Posted 2021-01-31 asuldb

这个题简直有毒，(O((a+b)^3logn))的做法不卡常只比(O(2^n*n))多(10)分

看到(a)和(b)简直小的可怜，于是可以往矩阵上联想

发现这个柿子有些特殊，好像可以二项式定理搞一搞

于是(x^ay^b)可以写成((n-y)^ay^b)

于是接下来就二项式定理好了

[(n-y)^ay^b=sum_{r=0}^ainom{a}{r}n^{a-r}*(-y)^r*y^b]

[=sum_{r=0}^ainom{a}{r}n^{a-r}*(-1)^r*y^{b+r}]

发现好像可以用矩阵来维护这个(sum)的每一项

先列一下(dp)的方程，设(dp[i][j][0/1])表示进行到第(i)位上，这个(sum)的第(j)次方项，最后一位填的是(0/1)

如果这一位填(0)，对答案并没有什么贡献，但是前面填(0/1)都是可以的，于是(dp[i][j][0]=dp[i-1][j][0]+dp[i-1][j][1])

如果这一位填的是(1)，那么前面的那一位只能填(0)，(y)增加了(1)，所以答案变成了((y+1)^j)

还是用二项式定理

[(y+1)^j=sum_{k=0}^jinom{j}{k}y^k]

所以也就可以得到

[dp[i][j][1]=sum_{k=0}^jinom{j}{k}*dp[i-1][k][0]]

矩阵维护就可以了

#include<iostream>
#include<cstring>
#include<cstdio>
#define re register
#define maxn 185
#define LL long long
LL a[maxn][maxn],ans[maxn][maxn];
int sz;
int T,A,b,P;
LL c[maxn][maxn];
inline void did_a()
{
    LL mid[maxn][maxn];
    for(re int i=1;i<=sz;i++)
        for(re int j=1;j<=sz;j++)
            mid[i][j]=a[i][j],a[i][j]=0;
    for(re int i=1;i<=sz;i++)
        for(re int k=1;k<=sz;k++)
            for(re int j=1;j<=sz;j++)
            {
                a[i][j]+=mid[i][k]*mid[k][j];
                if(a[i][j]>P) a[i][j]%=P;
            }
}
inline void did_ans()
{
    LL mid[maxn][maxn];
    for(re int i=1;i<=sz;i++)
        for(re int j=1;j<=sz;j++)
            mid[i][j]=ans[i][j],ans[i][j]=0;
    for(re int i=1;i<=sz;i++)
        for(re int k=1;k<=sz;k++)
            for(re int j=1;j<=sz;j++)
            {
                    ans[i][j]+=mid[i][k]*a[k][j];
                    if(ans[i][j]>P) ans[i][j]%=P;
            }
}
inline LL quick(LL a,int b)
{
    LL S=1;
    while(b)
    {
        if(b&1) S=S*a%P;
        b>>=1;
        a=a*a%P;
    }
    return S;
}
inline void Mat_quick(int b)
{
    while(b)
    {
        if(b&1) did_ans();
        b>>=1;
        did_a();
    }
}
int main()
{
    scanf("%d%d%d%d",&T,&A,&b,&P);
    c[0][0]=1;
    for(re int i=1;i<=A+b;i++) c[i][0]=c[i][i]=1;
    for(re int i=1;i<=A+b;i++)
        for(re int j=1;j<i;j++)
            c[i][j]=(c[i-1][j-1]+c[i-1][j])%P;
    sz=(A+b+1)*2;
    for(re int i=1;i<=sz;i++)
        ans[i][i]=1;
    for(re int i=1;i<=(sz>>1);i++)
        a[i][i]=a[i+sz/2][i]=1;
    for(re int i=1;i<=(sz>>1);i++)
        for(re int j=i+sz/2;j<=sz;j++)
            a[i][j]=c[j-1-sz/2][i-1];
    Mat_quick(T);
    LL Ans=0;
    for(re int i=0;i<=A;i++)
    if(i&1) Ans=(Ans-c[A][i]*quick(T,A-i)%P*(ans[1][i+1+b]+ans[1][i+1+b+sz/2]%P)%P+P)%P;
        else Ans=(Ans+c[A][i]*quick(T,A-i)%P*(ans[1][i+1+b]+ans[1][i+1+b+sz/2]%P))%P;
    std::cout<<Ans;
    return 0;
}