题解：洛谷P1357 花园
Description
小 L 有一座环形花园，沿花园的顺时针方向，他把各个花圃编号为 \(1∼n\)。花园 \(1\) 和 \(n\) 是相邻的。
他的环形花园每天都会换一个新花样，但他的花园都不外乎一个规则：任意相邻 \(m\) 个花圃中都只有不超过 \(k\) 个 C 形的花圃，其余花圃均为 P 形的花圃。
例如，若 \(n=10 , m=5 , k=3\) ，则
· CCPCPPPPCC 是一种不符合规则的花圃。
· CCPPPPCPCP 是一种符合规则的花圃。
请帮小 L 求出符合规则的花园种数对 \(1e9 + 7\) 取模的结果。
Hints
对于\(40\%\)的数据，\(n \leq 20\)
对于\(60\%\)的数据，\(m = 2\)
对于\(80\%\)的数据，\(n \leq 10 ^5\)
对于\(100\%\)的数据，\(2 \leq n \leq 10^{15},\ 2 \leq m \leq min(n,5),\ 1 \leq k < m\)
Algorithm
单看Description可能以为这是个纯数学题，然而上手就能发现这玩意不太可做。
其实，这题其实已经把算法写在脸上了。
首先可以观察到 \(m\) 的值很小，因此显然需要状态压缩。
接着容易发现 \(n\) 的量级对于状态压缩来说大得离谱，因此显然是矩阵加速递推做的。
不管计数DP有没有最优子结构，姑且也算DP吧，那么这题就是个矩阵加速的计数类环形状压DP。
接下来把递推式搞出来就完了。
先给出一些约定：
我们考虑把放 C 当做 1 ，放 P 当做 0 ，将相邻 \(m\) 位的摆花方法当成一个二进制数压缩。
因为转移时左移操作会带来麻烦，所以把二进数的位置与原串摆法相反。
所有位运算符号均采用 C++ 中的样式。
Step 1
先考虑普通的状压DP：
我们用 dp[i - 1][j] 表示考虑到第 \(i - 1\) 个位置，后 \(m\) 位顺次摆花的状态为 \(j\) 时的方案数。
如果把第 \(i\) 位设为 0 ，那么显然有新状态 j' = j >> 1 ，
如果把第 \(i\) 位设为 1 ，那么显然有新状态 j' = (j >> 1) | (1 << (m - 1)) ，
那么 dp[i][j] 显然就是dp[i - 1][j >> 1] 与 dp[i - 1][(j >> 1) | (1 << (m - 1))]求和得到的。
注意后者不一定合法，需要 __builtin_popcount() 判断一下。//这是个统计二进制数中 1 的个数的函数，奇快无比。
最后，其实这题环的特征其实并不明显，
对于前 \(m\) 位状态为 t 的情况来说，我们只需要递推完事之后将 dp[n][t] 计入答案即可。
当然，容易注意到 dp[i][] 只与 dp[i - 1][] 有关，可以滚动优化一下：
细节的小技巧就不细说了，容易写出代码：
int ans = 0, ful = 1 << m, p = m - 1;
for(int t = 0; t != ful; ++t)
{
if(__builtin_popcount(t) > k) continue;
memset(dp, 0, sizeof(dp));
dp[0][t] = 1;
for(int i = 1; i <= n; ++i)
for(int j = 0; j <= ful; ++j)
{
dp[i & 1][j] = dp[i - 1 & 1][j >> 1];
if(__builtin_popcount((j >> 1) | (1 << p)) <= k)
dp[i & 1][j] += dp[i - 1 & 1][(j >> 1) | (1 << p)];
}
ans += dp[n & 1][t];
}
复杂度显然是 \(O(n\cdot 4^m)\) 的。
Step 2
接下来就需要构造递推矩阵了。
以 \(m = 2, k = 1\) 为例，我们有
dp[n][00] = dp[n - 1][00] + dp[n - 1][01];
dp[n][01] = dp[n - 1][10];
dp[n][10] = dp[n - 1][00] + dp[n - 1][01];
dp[n][11] = 0;
从而有转移方程：
\begin{matrix}
dp_{n，00} \\
dp_{n，01} \\
dp_{n，10} \\
dp_{n，11} \\
\end{matrix}
\right]
=
\left[
\begin{matrix}
1 & 1 & 0 & 0 \\
0 & 0 & 1 & 0 \\
1 & 1 & 0 & 0 \\
0 & 0 & 0 & 0
\end{matrix}
\right]
\left[
\begin{matrix}
dp_{n - 1，00} \\
dp_{n - 1，01} \\
dp_{n - 1，10} \\
dp_{n - 1，11} \\
\end{matrix}
\right]
\]
（\(\LaTeX\) 的矩阵也太难写了吧……）
类似地，对于具体的 \(m,k\) ，
我们可以使用和上面的递推程序相同的方法构造出类似的矩阵。
最后同样来考虑环的问题。
根据 Step 1 中的分析，答案实际上就是\(\sum dp[\ n\ ][\ j\ ]\)，在我们这里又有
\begin{matrix}
dp_{n，00\dots0} \\
dp_{n，00\dots1} \\
\vdots \\
dp_{n，11\dots0} \\
dp_{n，11\dots1} \\
\end{matrix}
\right]
=
\begin{bmatrix}
1 & 1 & \cdots & 0 & 0 \\
0 & 0 & \cdots & 0 & 0 \\
\vdots & \vdots & \ddots & \vdots & \vdots \\
0 & 0 & \cdots\ & 0 & 0\\
0 & 0 & \cdots\ & 1 & 1 \\
\end{bmatrix}
^ n
\left[
\begin{matrix}
1\\
1\\
\vdots \\
0\\
0\\
\end{matrix}
\right]
\]
（特别应当注意，右侧列向量并非全 1 向量）
如果把右侧列向量表示出来乘上去再求和，自然能做，不过有些麻烦了。
这里实际上有个绝妙的性质：
考虑构造出来的矩阵 \(A\) , \(A_{i,j}\) 实际上表示了状态 \(j\) 可以向状态 \(i\) 转移。
而由于环首尾相连的特性，对于一个合法的状态（也就是右侧列向量中对应值为 1 的状态），
在 \(n\) 次转移后，该状态一定可以转移回原状态，也即是说 \(A_{i,i}\) 一定非零。
而在矩阵做了 \(n\) 次幂乘之后， \(A_{i,i}\) 也就是起始状态为 \(i\) 的方案数。
因此答案实际上即为 \(\sum A_{i,i}\)。
\(\sum A_{i,i} = \sum dp[\ n\ ][\ j\ ]\) 的代数证明留给读者。（这能证吗？不管了，我们读者肯定会）
最后写写代码：
#include<bits/stdc++.h>
using namespace std;
#define cnt __builtin_popcount
#define rep(i, n) for(int i = 0; i != n; ++i)
typedef long long ll;
const int P = 1e9 + 7;
ll n, m, k;
class Matrix {
private:
vector<vector<int> > num;
public:
Matrix() {}
Matrix(int n) {
num.resize(n);
rep(i, n) num[i].resize(n, 0);
}
int &operator () (int i, int j) {
return num[i][j];
}
inline void print()
{
int len = num.size();
for(int i = 0; i != len; ++i)
{
printf("%d", num[i][0]);
for(int j = 1; j != len; ++j)
printf(" %d", num[i][j]);
puts("");
}
}
inline int solve()
{
int ans = 0, len = num.size();
for(int i = 0; i != len; ++i)
(ans += num[i][i]) %= P;
return ans;
}
friend Matrix operator * (Matrix a, Matrix b)
{
int len = a.num.size();
Matrix ret(len);
rep(k, len) rep(i, len) rep(j, len) {
(ret(i, j) += (1LL * a(i, k) * b(k, j)) % P) %= P;
}
return ret;
}
friend Matrix operator ^ (Matrix a, ll n)
{
int len = a.num.size();
Matrix ret(len);
rep(i, len) ret(i, i) = 1;
while(n)
{
if(n & 1) ret = ret * a;
a = a * a, n >>= 1;
}
return ret;
}
};
int main()
{
cin >> n >> m >> k;
int ful = 1 << m, p = m - 1;
Matrix cal(ful);
for(int t = 0; t != ful; ++t)
{
if(cnt(t) > k) continue;
cal(t >> 1, t) = 1;
int sta = (t >> 1) | (1 << p);
if(cnt(sta) <= k) cal(sta, t) = 1;
}
cout << (cal ^ n).solve() << endl;
return 0;
}
这个代码非常之丑。最开始设计模板类翻车了，没有重构而是左改右改的……
总之不是好代码，不过反正是这个意思。
复杂度 \(O(log n \cdot 8 ^ m)\) 很容易过。