用 2 x 1 和 1 x 2 多米诺骨牌拼贴 W x H 网格的方法有多少种?
Number of ways to tile a W x H Grid with 2 x 1 and 1 x 2 dominos?
编写一个程序,将网格的宽度、W 和高度 H 作为输入,并输出用 (2x1) 个多米诺骨牌平铺 W×H 网格的不同方法的数量
我知道如何解决 3 x N 网格的问题,但是为此编写递归公式对我来说太难了!
我不知道该怎么做!
我创建了两个函数 F(n) - The complete way of tiling till N 和 S(n) for number of ways of incomplete tiling for 3 x N !但是这里因为高度是可变的我想不出任何东西
约束:(0≤W+H≤22)
有一种很好的通用方法可以解决这类问题:
1) 计算前几个值。
0 1 0 1 0 1
1 2 3 5 8 13
0 3 0 11 0 41
1 5 11 36 95 281
0 8 0 95 0 1183
1 13 41 281 1183 6728
2) 在整数序列在线百科全书中搜索这些值。您可以通过输入连接的反对角线来搜索方形数组,或者您可以对一些大项进行无序搜索。成功:https://oeis.org/A099390 其中包含涉及余弦的乘积公式。
据我所知,没有一种特别好的方法可以直接通过动态规划来执行此操作,比如在最小维度中使用少于指数级的状态。与一些高等数学有联系。请参阅 OEIS 条目或 Wikipedia.
中的参考资料
Kasteleyn 的余弦积公式可以用精确的算法来实现。其实可以用动态规划来完成。
2 cos x = e^ix + e^-ix
所以,Kasteleyn 的产品
(4+2 cos 2pi j/(w+1) + 2 cos 2pi k/(h+1))
在 1<= j <= w/2 范围内,1 <= k <= h/2 可以写成 z=e^i 2pi/((w+1)*(h+1))
中的多项式,其中可以计算指数mod (w+1)*(h+1)
因为 z^((w+1)*(h+1))=z^0。可以使用动态规划计算多项式的乘积,您可以跟踪从 0 到 (w+1)*(h+1)-1
的每个 n 的 z^n 系数。
// part of a C# implementation using a custom CircularVector class
using System.Numerics; // BigInteger
CircularVector cv = new CircularVector(w, h);
cv.setToOne();
for (int i = 1; i<=w/2; i++)
for (int j = 1; j <= h / 2; j++)
cv = cv.multiplyWith(new CircularVector(w, h, i, j));
结果是最多 (w+1)*(h+1)
项的总和,是 z 中小于 (w+1)*(h+1)
次的多项式。例如,w=3,h=4,我们得到
{18, 1, 0, 1, 5, 8, 0, 1, 5, 1, 2, 1, 5, 1, 0, 8, 5, 1, 0, 1}
= 18 + z + z^3 + 5z^4 + 8z^5 + z^7 + 5z^8 + z^9 + 2 z^10 + ...
需要找到等于此多项式 z 的整数。为此,我们可以将多项式除以 cyclotomic polynomial Phi((w+1)*(h+1))
.
// part of a C# implementation
public static BigInteger[] cyclotomic(int n)
{
BigInteger[] working = nthRoots(n); // x^n - 1
for (int i = 1; i<n; i++)
if (n % i == 0)
working = polynomialExactDivision(working, cyclotomic(i));
return working;
}
除以分圆多项式后的余数是(大)整数,计算矩形的多米诺骨牌拼贴数。这在一个处理器上花费了 0.014 秒来计算 20x20 正方形的平铺数量,1269984011256235834242602753102293934298576249856,这与 OEIS A004003 中的值一致。计算可以在 O(w^2h^2) BigInteger 加法中完成,因为您执行大小为 wh 的多项式乘以 5 项的因子的 wh 乘法。
这是引擎盖下的 DPish。这个想法是 IVlad 的:使用记忆回溯。
memo = {}
def count_tilings_recursive(uncovered):
if len(uncovered) & 1:
return 0
if not uncovered:
return 1
if uncovered not in memo:
i, j = min(uncovered)
memo[uncovered] = count_tilings_recursive(uncovered - {(i, j), (i, j + 1)}) + count_tilings_recursive(uncovered - {(i, j), (i + 1, j)})
return memo[uncovered]
def count_tilings(m, n):
return count_tilings_recursive(frozenset((i, j) for i in range(max(m, n)) for j in range(min(m, n))))
print(count_tilings(18, 4))
这里的诀窍是防止状态数量激增。首先,我们总是选择最左边然后最上面的方块进行覆盖,因此部分覆盖可以描述为按最左-最上面的顺序连续覆盖的方块的数量,然后是部分覆盖的#rows,总共最多(#rows * #columns + 1) * 2^#rows 状态。其次,我们调整网格的方向,使列数至少与行数一样多,将后者的数字限制为 10 以进行有趣的计算(因为 11 乘 11 是微不足道的)。
在代码中写出不同状态之间的关系。高度为2
的情况等价于斐波那契数列。高度 3
的情况可以用具有九个相互递归序列的公式来编写(有点乏味!),或者等效地使用大小为 9
.
的向量的线性递归
H>=4
,请看这个article。
R。 C. 阅读,关于用多米诺骨牌拼接矩形的注释,斐波那契季刊,18.1 (1980),24-27。
从n * 2
格子的情况说起。填充方式的数量是第n
-th斐波那契数。为什么?因为右端的格子有两种补全的方式,要么用一张垂直的多米诺骨牌,要么用两张水平的多米诺骨牌。
查看它的另一种方法是从空网格开始,然后从左到右填充它。考虑到这一点,我们可以根据它们的 "profile" 类型对所有可能的不完整网格进行分组,即右侧的边界。
例如,对于H=2
案例,有两种配置文件,对于H=3
,有九种配置文件
(文章中的图2)
对于一般情况,使用 H >= 2
,我们 pre-compute 配置文件列表。每个配置文件可以由两个位掩码编码。注意第一个不要代表一个完整的列,这样才能代表的没有歧义。
为了枚举配置文件类型,同时链接它们,我们尝试先水平放置一个多米诺骨牌,然后垂直放置在最左边、最低的可用单元格中。这会生成另一个配置文件,如果它是新的,我们将其添加到列表中。这样,每个配置文件最多可以链接到两个配置文件。
有关 H
的其他选择,请参阅下面的代码。这个过程为我们提供了一个递归 "formula" 在填充 r
多米诺骨牌的不完整网格和填充 r+1
多米诺骨牌的网格之间。这不是一个真正的公式,它被编码在 next_profile
字典中。相当于文章中的公式(1.1)
,对于H=3
.
from collections import defaultdict
def disambiguation(a,b):
if a == (1<<H) - 1:
return disambiguation(b,0)
return (a,b)
def generate_next_profile():
global next_profile
next_profile = {}
q = [(0,0)]
done = {(0,0): True}
while q != []:
a,b = q.pop(0)
next_profile[(a,b)] = []
i = 0
while i < H and a & 1<<i != 0:
i += 1
if i+1 < H and a & 1<<(i+1) == 0:
c, d = disambiguation(a | 1<<i | 1<<(i+1), b )
next_profile[(a,b)].append((c,d))
if (c,d) not in done:
q.append((c,d))
done[(c,d)] = True
c, d = disambiguation(a | (1<<i), b | (1<<i))
next_profile[(a,b)].append((c,d))
if (c,d) not in done:
q.append((c,d))
done[(c,d)] = True
def count_it(_W, _H):
W = _W
H = _H
result = 0
if W * H % 2 == 0:
u = 0
v = 1
num = [{(0,0): 1}, {}]
n_steps = W * H // 2
for step in range(n_steps):
num[v] = defaultdict(int)
for key, value in num[u].items():
for pair in next_profile[key]:
num[v][pair] += value
u, v = v, u
result = num[n_steps % 2][(0,0)]
return result
WMAX = 10
HMAX = 10
for W in range(1, WMAX):
for H in range(1, HMAX):
generate_next_profile()
res = count_it(W, H)
print(res, end = "\t")
print()
编写一个程序,将网格的宽度、W 和高度 H 作为输入,并输出用 (2x1) 个多米诺骨牌平铺 W×H 网格的不同方法的数量
我知道如何解决 3 x N 网格的问题,但是为此编写递归公式对我来说太难了!
我不知道该怎么做!
我创建了两个函数 F(n) - The complete way of tiling till N 和 S(n) for number of ways of incomplete tiling for 3 x N !但是这里因为高度是可变的我想不出任何东西
约束:(0≤W+H≤22)
有一种很好的通用方法可以解决这类问题:
1) 计算前几个值。
0 1 0 1 0 1
1 2 3 5 8 13
0 3 0 11 0 41
1 5 11 36 95 281
0 8 0 95 0 1183
1 13 41 281 1183 6728
2) 在整数序列在线百科全书中搜索这些值。您可以通过输入连接的反对角线来搜索方形数组,或者您可以对一些大项进行无序搜索。成功:https://oeis.org/A099390 其中包含涉及余弦的乘积公式。
据我所知,没有一种特别好的方法可以直接通过动态规划来执行此操作,比如在最小维度中使用少于指数级的状态。与一些高等数学有联系。请参阅 OEIS 条目或 Wikipedia.
中的参考资料Kasteleyn 的余弦积公式可以用精确的算法来实现。其实可以用动态规划来完成。
2 cos x = e^ix + e^-ix
所以,Kasteleyn 的产品
(4+2 cos 2pi j/(w+1) + 2 cos 2pi k/(h+1))
在 1<= j <= w/2 范围内,1 <= k <= h/2 可以写成 z=e^i 2pi/((w+1)*(h+1))
中的多项式,其中可以计算指数mod (w+1)*(h+1)
因为 z^((w+1)*(h+1))=z^0。可以使用动态规划计算多项式的乘积,您可以跟踪从 0 到 (w+1)*(h+1)-1
的每个 n 的 z^n 系数。
// part of a C# implementation using a custom CircularVector class
using System.Numerics; // BigInteger
CircularVector cv = new CircularVector(w, h);
cv.setToOne();
for (int i = 1; i<=w/2; i++)
for (int j = 1; j <= h / 2; j++)
cv = cv.multiplyWith(new CircularVector(w, h, i, j));
结果是最多 (w+1)*(h+1)
项的总和,是 z 中小于 (w+1)*(h+1)
次的多项式。例如,w=3,h=4,我们得到
{18, 1, 0, 1, 5, 8, 0, 1, 5, 1, 2, 1, 5, 1, 0, 8, 5, 1, 0, 1}
= 18 + z + z^3 + 5z^4 + 8z^5 + z^7 + 5z^8 + z^9 + 2 z^10 + ...
需要找到等于此多项式 z 的整数。为此,我们可以将多项式除以 cyclotomic polynomial Phi((w+1)*(h+1))
.
// part of a C# implementation
public static BigInteger[] cyclotomic(int n)
{
BigInteger[] working = nthRoots(n); // x^n - 1
for (int i = 1; i<n; i++)
if (n % i == 0)
working = polynomialExactDivision(working, cyclotomic(i));
return working;
}
除以分圆多项式后的余数是(大)整数,计算矩形的多米诺骨牌拼贴数。这在一个处理器上花费了 0.014 秒来计算 20x20 正方形的平铺数量,1269984011256235834242602753102293934298576249856,这与 OEIS A004003 中的值一致。计算可以在 O(w^2h^2) BigInteger 加法中完成,因为您执行大小为 wh 的多项式乘以 5 项的因子的 wh 乘法。
这是引擎盖下的 DPish。这个想法是 IVlad 的:使用记忆回溯。
memo = {}
def count_tilings_recursive(uncovered):
if len(uncovered) & 1:
return 0
if not uncovered:
return 1
if uncovered not in memo:
i, j = min(uncovered)
memo[uncovered] = count_tilings_recursive(uncovered - {(i, j), (i, j + 1)}) + count_tilings_recursive(uncovered - {(i, j), (i + 1, j)})
return memo[uncovered]
def count_tilings(m, n):
return count_tilings_recursive(frozenset((i, j) for i in range(max(m, n)) for j in range(min(m, n))))
print(count_tilings(18, 4))
这里的诀窍是防止状态数量激增。首先,我们总是选择最左边然后最上面的方块进行覆盖,因此部分覆盖可以描述为按最左-最上面的顺序连续覆盖的方块的数量,然后是部分覆盖的#rows,总共最多(#rows * #columns + 1) * 2^#rows 状态。其次,我们调整网格的方向,使列数至少与行数一样多,将后者的数字限制为 10 以进行有趣的计算(因为 11 乘 11 是微不足道的)。
在代码中写出不同状态之间的关系。高度为2
的情况等价于斐波那契数列。高度 3
的情况可以用具有九个相互递归序列的公式来编写(有点乏味!),或者等效地使用大小为 9
.
H>=4
,请看这个article。
R。 C. 阅读,关于用多米诺骨牌拼接矩形的注释,斐波那契季刊,18.1 (1980),24-27。
从n * 2
格子的情况说起。填充方式的数量是第n
-th斐波那契数。为什么?因为右端的格子有两种补全的方式,要么用一张垂直的多米诺骨牌,要么用两张水平的多米诺骨牌。
查看它的另一种方法是从空网格开始,然后从左到右填充它。考虑到这一点,我们可以根据它们的 "profile" 类型对所有可能的不完整网格进行分组,即右侧的边界。
例如,对于H=2
案例,有两种配置文件,对于H=3
,有九种配置文件
(文章中的图2)
对于一般情况,使用 H >= 2
,我们 pre-compute 配置文件列表。每个配置文件可以由两个位掩码编码。注意第一个不要代表一个完整的列,这样才能代表的没有歧义。
为了枚举配置文件类型,同时链接它们,我们尝试先水平放置一个多米诺骨牌,然后垂直放置在最左边、最低的可用单元格中。这会生成另一个配置文件,如果它是新的,我们将其添加到列表中。这样,每个配置文件最多可以链接到两个配置文件。
有关 H
的其他选择,请参阅下面的代码。这个过程为我们提供了一个递归 "formula" 在填充 r
多米诺骨牌的不完整网格和填充 r+1
多米诺骨牌的网格之间。这不是一个真正的公式,它被编码在 next_profile
字典中。相当于文章中的公式(1.1)
,对于H=3
.
from collections import defaultdict
def disambiguation(a,b):
if a == (1<<H) - 1:
return disambiguation(b,0)
return (a,b)
def generate_next_profile():
global next_profile
next_profile = {}
q = [(0,0)]
done = {(0,0): True}
while q != []:
a,b = q.pop(0)
next_profile[(a,b)] = []
i = 0
while i < H and a & 1<<i != 0:
i += 1
if i+1 < H and a & 1<<(i+1) == 0:
c, d = disambiguation(a | 1<<i | 1<<(i+1), b )
next_profile[(a,b)].append((c,d))
if (c,d) not in done:
q.append((c,d))
done[(c,d)] = True
c, d = disambiguation(a | (1<<i), b | (1<<i))
next_profile[(a,b)].append((c,d))
if (c,d) not in done:
q.append((c,d))
done[(c,d)] = True
def count_it(_W, _H):
W = _W
H = _H
result = 0
if W * H % 2 == 0:
u = 0
v = 1
num = [{(0,0): 1}, {}]
n_steps = W * H // 2
for step in range(n_steps):
num[v] = defaultdict(int)
for key, value in num[u].items():
for pair in next_profile[key]:
num[v][pair] += value
u, v = v, u
result = num[n_steps % 2][(0,0)]
return result
WMAX = 10
HMAX = 10
for W in range(1, WMAX):
for H in range(1, HMAX):
generate_next_profile()
res = count_it(W, H)
print(res, end = "\t")
print()