「题解」摆花

矩阵快速幂优化二维DP……

一 题目

原题链接

描述

艺术馆门前将摆出许多花，一共有n个位置排成一排，每个位置可以摆花也可以不摆花。有些花如果摆在相邻的位置（隔着一个空的位置不算相邻），就不好看了。假定每种花数量无限，求摆花的方案数。

输入

输入有1+m行，第一行有两个用空格隔开的正整数n、m，m表示花的种类数。接下来的m行，每行有m个字符1或0,若第i行第j列为1，则表示第i种花和第j种花不能排在相邻的位置，输入保证对称。（提示：同一种花可能不能排在相邻位置）。

输出

输出只有一个整数，为方案数（这个数字可能很大，请输出方案数除以1000000007的余数）。

---

二 题解

对于矩阵的初学者，建议您先看矩阵基础。

DP思路

很容易想到一个$O(nm)$级别的DP：
设$f_{i,j}$为考虑前$i$个位置，且第$i$个位置摆放第$j$种花时的方案数。由于任意位置都可以一盆花也不摆，故当$j = 0 $时表示该位置不摆花。

初始状态： $f_{0,0} = 1$
转移方程：设$a(x,y) = 1$表示$x,y$两种花可以相邻。则转移方程为： $ f_{i,j} = \sum\limits_{k\in [0,m]}^{a(i,k) = 1} f_{i-1,k}$
答案： 答案为$\sum\limits_{i = 0}^m f_{n,i}$

空间复杂度为$O(nm)$，朴素时间复杂度为$O(nm^2)$。核心代码大约长这样：

const int R = 1e9+7;
f[0][0]=1;
for(int i=1;i<=n;i++)
    for(int j=0;j<=m;j++)
        for(int k=0;k<=m;k++)
            if(a[j][k] || j==0 || k==0)
                f[i][j] = (f[i][j]+f[i-1][k])%R;
for(int i=0;i<=m;i++) ans = (ans+f[n][i])%R;

矩阵加速

我们发现：上面的这个DP转移是由前一层状态转移到后一层状态，而我们最终需要第$n$层状态。

我们不妨设$F(i)$表示第$i$层的状态，即$F(i) = \begin{bmatrix} f_{i,0}&f_{i,1}&f_{i,2}&…&f_{i,m} \end{bmatrix}$，这是一个$1\times (m+1)$阶的矩阵。

显然有$F(0) = \begin{bmatrix}1&\overbrace{0 \text{ }\text{ }\text{ }\text{ }0\text{ }\text{ }\text{ }\text{ }…} \end{bmatrix}$，后面的大括号里一共有$m$个$0$。

现在我们考虑怎么由$F(i-1)$推导$F(i)$。不妨把$a(x,y)$函数看成一个$(m+1)\times (m+1)$阶的矩阵$A$，则有（注：乘号后面的$ A_{0\text{~}m,j}$表示$A$矩阵的第$j$列，即一个$1\times (m+1)$阶的矩阵）：
$$ \begin{aligned}& F(i-1) \times A_{0\text{~}m,0} = f_{i,0} \newline & F(i-1) \times A_{0\text{~}m,1} = f_{i,1} \newline &… \newline & F(i-1) \times A_{0\text{~}m,m} = f_{i,m} \end{aligned}$$

故：$ F(i-1) \times A = F(i)$$ $$ F(i) \times A^k = F(i+k)$

那么我们就可以用$F(0) \times A^n = F(n)$推导出第$n$层状态，再求和就好了。

时间复杂度大约是$O(m^2\log n )$。

一些细节
$0$可以和任意一种花搭配，故$F,A$矩阵都是$k\times {m+1}$阶级别的。因此在代码里我把$A$矩阵强行向右下平移了一格，即矩阵的有效区间是$(1\text{~}m+1) \times (1 \text{~}m+1) $。在代码里看可能会更直白一些。

又因为$0$可以和任意一种花搭配，所以$A_{1,i} = A_{i,1} = 1 (i\in[1,m+1])$（矩阵已经平移过了）。一定注意这一点！

代码：

#include<iostream>
#include<cstdio>
#include<cstring>
#include<string>
using namespace std;

#define LL long long

const int CM = 102;
const LL R = 1e9+7;

int n,m;

class matrix{ //定义矩阵
  public: LL a[CM][CM]; int n,m;
    void MakeOne(int k){
        n = m = k;
        for(int i=1;i<=n;i++) a[i][i] = 1;
    }
    matrix operator * (const matrix &b)const{
        matrix rec = (matrix){{{0}}, n, b.m}; //注意别忘了定义相乘得到的矩阵的大小
        for(int i=1;i<=n;i++)
            for(int j=1;j<=m;j++)
                for(int k=1;k<=b.m;k++)
                    (rec.a[i][k] += (a[i][j]*b.a[j][k])%R ) %= R; 
        return rec;
    }
    void operator *= (const matrix &b) {*this = (*this)*b;}
}A;

matrix MatrixQuickPow(matrix &a,LL b){ //矩阵快速幂
    matrix rec; rec.MakeOne(a.m);
    while(b){
        if(b & 1) rec *= a;
        a *= a; b >>= 1;
    }
    return a = rec;
}

int main()
{
    //freopen("in.in","r",stdin);
    
    scanf("%d%d",&n,&m);
    A.n = A.m = m+1; //平移后矩阵的大小是 (m+1)^2
    
    for(int i=1;i<=m;i++)
        for(int j=1;j<=m;j++){
            char c; cin>>c;
            if(c == '0') A.a[i+1][j+1] = 1; //平移
        }
    for(int i=1;i<=m+1;i++) A.a[1][i] = A.a[i][1] = 1; //初始化，0可以和任意花搭配
        
    matrix c = (matrix){{{0}}, 1, m+1}; //定义 F(0) 
    c.a[1][1] = 1; //初始状态
    c *= MatrixQuickPow(A, n); //递推到 F(n)
    
    LL ans = 0;
    for(int i=1;i<=m+1;i++) (ans += c.a[1][i]) %= R; //累加答案
    
    printf("%lld",ans);
    
    return 0;
}