多米诺与托米诺骨牌平铺:矩阵快速幂优化DP

  |  

摘要: 动态规划解决计数问题,模下矩阵快速幂优化

【对数据分析、人工智能、金融科技、风控服务感兴趣的同学,欢迎关注我哈,阅读更多原创文章】
我的网站:潮汐朝夕的生活实验室
我的公众号:潮汐朝夕
我的知乎:潮汐朝夕
我的github:FennelDumplings
我的leetcode:FennelDumplings

各位好,之前在文章 计数DP 中我们总结了常见的使用动态规划解决计数问题的题目,在文章 网格图涂色的方案数 中我们解决了一个使用 3 种颜色给 3 * n 网格图涂色方案数的问题。

本文我们来看一个用动态规划解决计数问题的例子:骨牌平铺问题。本题的状态转移方程式线性递推关系,可以写成矩阵形式,然后用矩阵快速幂优化,最后的代码可以作为模下矩阵快速幂的代码模板使用。

题目

790. 多米诺和托米诺平铺

有两种形状的瓷砖:一种是 2 x 1 的多米诺形,另一种是形如 “L” 的托米诺形。两种形状都可以旋转。

给定整数 n ,返回可以平铺 2 x n 的面板的方法的数量。返回对 1e9 + 7 取模 的值。

平铺指的是每个正方形都必须有瓷砖覆盖。两个平铺不同,当且仅当面板上有四个方向上的相邻单元中的两个,使得恰好有一个平铺有一个瓷砖占据两个正方形。

提示:

1
1 <= n <= 1000

示例 1:
输入: n = 3
输出: 5
解释: 五种不同的方法如上所示。

示例 2:
输入: n = 1
输出: 1

题解

算法1:动态规划

定义 f[i] 为铺满 2 * i 的面板的方法数。f[n] 为我们要求的值。

首先看 i = 1, 2 的情况:f[1] = 1, f[2] = 2。这两种情况比较好想,列数为 1 时,只能使用多米诺型骨牌,并且只能竖着放;当列数为 2 时,也是只能使用多米诺型骨牌,但是可以两块都横着放或者两块都竖着放。

当 i >= 3 的时候,就可以使用托米诺型骨牌了。考虑最右边的一列怎么放:

可以竖着放多米诺骨牌,则剩下 (i - 1) * 2 的面板待铺满;

也可以横着放多米诺骨牌,这样需要占用两列,则剩下 (i - 2) * 2 的面板待铺满;

还可以用托米诺型骨牌,这会使得剩下的部分会在 (i - 2) * 2 的面板的基础上,往外凸出一个位置。为了解决这个情况,我们额外定义 g[i]i * 2 面板基础上向外凸出一块时,铺满的方法数。

同样地,首先写出 i = 1, 2 的情况:g[1] = 1, g[2] = 2,入上图。

定义好 g[i] 后,我们就可以完整推导 f[i] 的转移方程了,考虑最右边一列的放法,可以分成四种情况,如下图:

得到 f[i] 状态转移方程:f[i] = f[i - 1] + f[i - 2] + 2 * g[i - 2]。然后我推导 g[i] 的状态转移方程,如下图:

得到 g[i] 的状态转移方程:g[i] = g[i - 1] + f[i - 1]

综上,写出动态规划算法如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
状态定义:
f[i] := 铺满 `2 * i` 的面板的方法数
g[i] := `i * 2` 面板基础上向外凸出一块时,铺满的方法数。

答案:
f[n]

初始化:
f[1] = 1, f[2] = 2
g[1] = 1, g[2] = 2

状态转移:
f[i] = f[i - 1] + f[i - 2] + 2 * g[i - 2]
g[i] = g[i - 1] + f[i - 1]

代码 (C++)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
class Solution {
public:
int numTilings(int N) {
if(N == 1) return 1;
if(N == 2) return 2;
vector<int> f(N + 1, 0), g(N + 1, 0);
f[1] = 1;
f[2] = 2;
g[1] = 1;
g[2] = 2;
for(int i = 3; i <= N; ++i)
{
f[i] = ((ll)f[i - 1] + (ll)f[i - 2] + 2 * (ll)g[i - 2]) % MOD;
g[i] = ((ll)g[i - 1] + (ll)f[i - 1]) % MOD;
}
return f[N];
}

private:
using ll = long long;
const int MOD = 1e9 + 7;
};

优化:矩阵快速幂

前面的动态规划算法的状态转移方程如下:

由于这个状态转移方程式线性递推的关系,因此我们可以考虑将这种线性递推关系用矩阵的形式写出来。

首先做一些推导,把 $g(i)$ 去掉。

由第一个式子,写出以下两个式子,分别记为【式子 1】和【式子 2】:

由第二个式子,写出下面一个式子,记为【式子 3】:

用【式子 3】加上【式子 3】再减去【式子 3】,得到以下公式:

也就是:

将这个式子写成用矩阵表示的矩阵的递推形式:

根据前面的线性递推关系,将 A 中的元素填充如下:

事先计算好 f[1..4]: $f(1) = 1$, $f(2) = 2$, $f(3) = 5$, $f(4) = 11$。

于是最终结果如下:

代码 (C++) 模下矩阵快速幂的代码模板

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
using ll = long long;
const int MOD = 1e9 + 7;

class Matrix {
public:
Matrix(int m){
M = m;
mat = vector<vector<int> >(M, vector<int>(M, 0));
}

void assign(int i, int j, int val)
{
mat[i][j] = val;
}

int get(int i, int j)
{
return mat[i][j];
}

void init()
{
for(int i = 0; i < M; ++i)
for(int j = 0; j < M; ++j)
mat[i][j] = 0;
for(int i = 0; i < M; ++i)
mat[i][i] = 1;
}

Matrix operator*(const Matrix& matrix) const
{
Matrix ans(M);
for(int i = 0; i < M; ++i)
for(int j = 0; j < M; ++j)
for(int k = 0; k < M; ++k)
ans.mat[i][j] = ((ll)ans.mat[i][j] + ((ll)mat[i][k] * matrix.mat[k][j]) % MOD) % MOD;
return ans;
}

Matrix operator^(int n) const
{
Matrix ans(M);
ans.init();
Matrix A = *this; // 拷贝出来用于自乘
while(n)
{
if(n & 1)
ans = ans * A;
A = A * A;
n >>= 1;
}
return ans;
}

void operator=(const Matrix& matrix)
{
for(int i = 0; i < M; ++i)
for(int j = 0; j < M; ++j)
mat[i][j] = matrix.mat[i][j];
}

void show()
{
cout << "matrix" << endl;
for(int i = 0; i < M; ++i)
{
for(int j = 0; j < M; ++j)
cout << a[i][j] << " ";
cout << endl;
}
}

private:
vector<vector<int>> mat;
int M;
};

class Solution {
public:
int numTilings(int N) {
// f(i) = 2 f(i - 1) + f(i - 3)
if(N == 1) return 1;
if(N == 2) return 2;
if(N == 3) return 5;
if(N == 4) return 11;
Matrix mat(4);
mat.assign(0, 0, 2);
mat.assign(0, 2, 1);
mat.assign(1, 0, 1);
mat.assign(2, 1, 1);
mat.assign(3, 2, 1);
Matrix mat_power = mat ^ (N - 4);
int ans = ((ll)mat_power.get(0, 0) * 11 % MOD) % MOD;
ans = (ans + (ll)mat_power.get(0, 1) * 5 % MOD) % MOD;
ans = (ans + (ll)mat_power.get(0, 2) * 2 % MOD) % MOD;
ans = (ans + (ll)mat_power.get(0, 3) % MOD) % MOD;
return ans;
}
};

Share