Cf1151f dp+矩阵快速幂

0有x个，1有y个；定义dp[i][j]:i次时前x位置有j个1的概率；可以从加减不变转移；k很大用矩阵快速幂

#include <iostream>
#include <cstdio>
#include <cstring>
#include <string>
#include <algorithm>
#include <cmath>

using namespace std;

long long n,k;
int a[105];
int N=100;

long long mod=1e9+7;

struct node
{
    long long mat[105][105];
    friend node operator*(node x,node y)
    {
        node c;
        for(int i=0;i<=N;i++)
        {
            for(int j=0;j<=N;j++)
            {
                c.mat[i][j]=0;
                for(int k=0;k<=N;k++)
                {
                    c.mat[i][j]+=(x.mat[i][k]*y.mat[k][j])%mod;
                }
                c.mat[i][j]%=mod;
            }
        }
        return c;
    }
}sta,p,anss;

node powd(node x,long long y)
{
    node ans;
    for(int i=0;i<=N;i++)
        for(int j=0;j<=N;j++)
            if(i==j)
                ans.mat[i][j]=1;
            else
                ans.mat[i][j]=0;

    while(y>0)
    {
        if(y%2==1)
            ans=ans*x;
        x=x*x;
        y/=2;
    }
    return ans;
}

long long powd(long long x,long long y)
{
    long long ans=1;
    while(y>0)
    {
        if(y%2==1)
            ans=ans*x%mod;
        x=x*x%mod;
        y/=2;
    }
    return ans;
}
long long afac(long long x)
{
    return powd(x,mod-2);
}

int main() {

    while(~scanf("%lld%lld",&n,&k))
    {
        int sum=0;
        for(int i=1;i<=n;i++)
            scanf("%d",&a[i]),sum+=a[i];

        long long x=n-sum,y=sum;
        long long sum1=0,sum2=0;
        for(int i=1;i<=x;i++)
        {
            if(a[i]==1)
                sum1++;
        }
        sum2=sum-sum1;

        long long ac=afac(n*(n-1)/2);

        for(int i=0;i<=N;i++)
        {
            for(int j=0;j<=N;j++)
            {
                if(j==0 && i==sum1)
                    sta.mat[i][j]=1;
                else
                    sta.mat[i][j]=0;

                p.mat[i][j]=0;
                if(i>x || j>x)
                    continue;

                long long p1,p2;
                p1=(x-j)*(sum-j)%mod;p1=p1*ac%mod;
                p2=j*(j)%mod;p2=p2*ac%mod;

                if(j==i)
                    p.mat[i][j]=(1-p1-p2)%mod;
                else if(j==i-1)
                    p.mat[i][j]=p1;
                else if(j==i+1)
                    p.mat[i][j]=p2;
            }
        }

        anss=powd(p,k)*sta;
        long long tans=anss.mat[0][0]%mod;
        tans=(tans+mod)%mod;
        printf("%lld\n",tans);
    }
    return 0;
}