题解:洛谷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;

从而有转移方程:

\[\left[
\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\ ]\),在我们这里又有

\[\left[
\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)\) 很容易过。

最新文章

  1. sublime一些快捷键
  2. 谈谈php里的IOC控制反转,DI依赖注入
  3. 四则运算appNABCD模型
  4. Redis 连接问题
  5. ThreadLocal原理与模拟
  6. 为Go Web App 创建一个主页面
  7. STL整理
  8. SQL 根据指定字符拆分字符串
  9. mockjs学习总结(方便前端模拟数据,加快开发效率)
  10. 使用 EF Power Tool Code Frist 生成 Mysql 实体
  11. 再硬写一个最简单的HTTPSERVER
  12. 关于javascript 数组的正态分布排序的一道面试题
  13. 《C语言及程序设计初步》网络课程主页
  14. C++对C的函数拓展 - 占位参数
  15. JAVA集合类兄妹List和Set
  16. UVA 10618 Tango Tango Insurrection
  17. 【ApplicationContext】通过实现ApplicationContextAware接口获取bean
  18. OAuth 2.0 Salesforce &amp; Azure
  19. python2和python3网络访问包
  20. python-迭代器与生成器的区别

热门文章

  1. axure rp extension for chrome怎么用
  2. JAVA读取文件夹大小
  3. Appium之启动第一个App
  4. nmap端口扫描工具下载和安装使用
  5. pycharm写的代码提交到git上,提示需要merge失败时解决办法
  6. linux定时删除过期文件
  7. yum管理——linux字符界面安装图形化及两种界面的切换(3)
  8. python的多行注释
  9. Redis散列(Hash)的相关命令
  10. session深入探讨