Catalan 数之Python演示
关于Catalan 数,英文的下面网页可以参考:
https://brilliant.org/wiki/catalan-numbers/
以下卡特兰数(Catalan number)简介来自下面网页,Python演示代码为原创。
http://blog.miskcoo.com/2015/07/catalan-number
卡特兰数(Catalan number)是组合数学中一个重要的计数数列,在很多看似毫不相关地方都能见到它的身影
带限制条件的路径总数
首先我们来看一个问题:
在一个平面直角坐标系中,只能往右或往上走一个单位长度,问有多少种不同的路径可以从左下角
(
1
,
1
)
(1, 1)
(1,1) 走到右上角
(
n
,
n
)
(n, n)
(n,n),并且要求路径不能经过直线
y
=
x
y = x
y=x 上方的点,下图中的路径都是合法的(图片来源 Wikipedia)
如果没有限制条件,那么从左下角走到右上角一共有
2
n
2n
2n 步,有
n
n
n 步是往右,另外
n
n
n 步是往上,那么路径方案数就是
2
n
2n
2n 步中选择
n
n
n步往右,一共有
(
2
n
n
)
{2n} \choose {n}
(n2n) (即
C
2
n
n
C_{2n}^n
C2nn)种方案。
那么我们考虑一下这里面有多少种方案是不合法的。
首先对于每一种不合法的方案,它的路径一定与 y = x + 1 y = x + 1 y=x+1 有交。我们找到它与 y = x + 1 y = x + 1 y=x+1 的第一个交点,然后将这个点后面部分的路径关于 y = x + 1 y = x + 1 y=x+1 做一个对称。由于原来路径到达 ( n , n ) (n, n) (n,n),新的对称之后的路径就会到达 ( n − 1 , n + 1 ) (n - 1, n +1) (n−1,n+1)。这样我们把每一种不合法方案都对应到了一条从 ( 1 , 1 ) (1, 1) (1,1) 到 ( n − 1 , n + 1 ) (n - 1, n + 1) (n−1,n+1)的路径。
现在再来看是否每一条这样的路径都能对应到一种不合法方案,如果是,那么这就建立了一个一一映射的关系,也就是它们的方案总数相同。这是肯定的,因为每一条这样的路径必定与 y = x + 1 y = x + 1 y=x+1 有交,那么对称回去,就得到一条不合法方案。
由于从
(
1
,
1
)
(1, 1)
(1,1) 到
(
n
−
1
,
n
+
1
)
(n - 1, n + 1)
(n−1,n+1) 的路径有
(
n
−
1
+
n
+
1
n
−
1
)
{n - 1 + n + 1} \choose {n - 1}
(n−1n−1+n+1)
条,那么合法的方案就是
C n = ( 2 n n ) − ( 2 n n − 1 ) = 1 n + 1 ( 2 n n ) C_n = { {2n} \choose {n} } - { {2n} \choose {n - 1} } = \frac{1}{n + 1} { {2n} \choose {n} } Cn=(n2n)−(n−12n)=n+11(n2n)
我们把这个方案数记为 C n C_n Cn,这就是著名的 Catalan 数
我们来看看它的前几项( n n n 从 0 0 0 开始)
1, 1, 2, 5, 14, 42, 132, 429, 1430, 4862, 16796, 58786, 208012, 742900,
2674440, 9694845, 35357670, 129644790, 477638700, 1767263190
Python演示代码说明
合法路径要求不越过对角线的本质意义,就是要求合法路径的每一步的坐标 x ≥ y x≥y x≥y,也即向右的步数必须大于或等于向上的步数。类比到计算机动态内存管理中,就是进栈动作数量必须大于或等于出栈动作数量,这是必然的。在括号配对问题中,则是左括号必须必须大于或等于右括号的数量,lisp语言中使用一系列的括号来区分代码的层次结构,非常能体现Catalan 数的意义。语言结构层次可以通过括号来区分,也可以通过二叉树来展示,括号配对与二叉树是可以互相映射的。
为了简单,代码演示中,选择了最简单的“田”字格,即 n = 2 n=2 n=2,从(0,0)出发,只能向左或向上,一共有 ( 2 ∗ 2 2 ) = 6 {{2*2} \choose {2}}=6 (22∗2)=6种路径方案到达(2,2),但仅有2种方案未越过对角线 x = y x=y x=y,属于合法路径。这不难理解。
Catalan 数理解中有难度的地方,就是如何计算非法路径的数量。如上面转引文字所述,就是把每一条非法路径的越线部分,进行翻转,使终点到达
(
n
−
1
,
n
+
1
)
(n - 1, n +1)
(n−1,n+1),就是走“目”字格。目字格,显然只有4条路径,我们可以通过的横线从下到上编号
1
、
2
、
3
、
4
1、2、3、4
1、2、3、4。数学家们巧妙的地方就是发现目字格的4条路径与田字格的4条非法路径是一一对应的。
演示图中,天蓝色的斜线是翻转对称轴 y = x + 1 y = x + 1 y=x+1。普通路径画为绿色,需要翻转的部分画为黄色,对称翻转后画为紫色。
第1条、第2条是合法路径,
[
[
0
,
0
]
,
[
1
,
0
]
,
[
2
,
0
]
,
[
2
,
1
]
,
[
2
,
2
]
]
[[0, 0], [1, 0], [2, 0], [2, 1], [2, 2]]
[[0,0],[1,0],[2,0],[2,1],[2,2]]
[
[
0
,
0
]
,
[
1
,
0
]
,
[
1
,
1
]
,
[
2
,
1
]
,
[
2
,
2
]
]
[[0, 0], [1, 0], [1, 1], [2, 1], [2, 2]]
[[0,0],[1,0],[1,1],[2,1],[2,2]]
每一步都是横坐标大于纵坐标,没有越过对角线。
第3条是非法路径,
[
[
0
,
0
]
,
[
1
,
0
]
,
[
1
,
1
]
,
[
1
,
2
]
,
[
2
,
2
]
]
[[0, 0], [1, 0], [1, 1], [1, 2], [2, 2]]
[[0,0],[1,0],[1,1],[1,2],[2,2]]
因为第3步,就是第4个坐标点,横坐标小于纵坐标,
1
<
2
1<2
1<2。从该点以后的路径需要以
y
=
x
+
1
y = x + 1
y=x+1为对称轴,进行翻转,以映射到目字格第1条路径。
第4条是非法路径,
[
[
0
,
0
]
,
[
0
,
1
]
,
[
1
,
1
]
,
[
2
,
1
]
,
[
2
,
2
]
]
[[0, 0], [0, 1], [1, 1], [2, 1], [2, 2]]
[[0,0],[0,1],[1,1],[2,1],[2,2]]
因为第1步,就是第2个坐标点,横坐标小于纵坐标,
0
<
1
0<1
0<1。从该点以后的路径需要以
y
=
x
+
1
y = x + 1
y=x+1为对称轴,进行翻转,以映射到目字格第4条路径。
第5条是非法路径,
[
[
0
,
0
]
,
[
0
,
1
]
,
[
1
,
1
]
,
[
1
,
2
]
,
[
2
,
2
]
]
[[0, 0], [0, 1], [1, 1], [1, 2], [2, 2]]
[[0,0],[0,1],[1,1],[1,2],[2,2]]
因为第1步,就是第2个坐标点,横坐标小于纵坐标,
0
<
1
0<1
0<1。从该点以后的路径需要以
y
=
x
+
1
y = x + 1
y=x+1为对称轴,进行翻转,以映射到目字格第3条路径。
第6条是非法路径,
[
[
0
,
0
]
,
[
0
,
1
]
,
[
0
,
2
]
,
[
1
,
2
]
,
[
2
,
2
]
]
[[0, 0], [0, 1], [0, 2], [1, 2], [2, 2]]
[[0,0],[0,1],[0,2],[1,2],[2,2]]
因为第1步,就是第2个坐标点,横坐标小于纵坐标,
0
<
1
0<1
0<1。从该点以后的路径需要以
y
=
x
+
1
y = x + 1
y=x+1为对称轴,进行翻转,以映射到目字格第2条路径。
求一点关于直线的对称点公式见下面网页,
https://tieba.baidu.com/p/5198944956
若点 A 1 ( x 1 , y 1 ) A_1(x_1,y_1) A1(x1,y1)关于直线 y = k x + b y=kx+b y=kx+b成轴对称,则对应点 A 2 ( x 2 , y 2 ) A_2(x_2,y_2) A2(x2,y2)的坐标为:
x 2 = 2 ⋅ k ⋅ ( y 1 − b ) − x 1 ⋅ ( k 2 − 1 ) k 2 + 1 x_2=\frac{2 \cdot k \cdot\left(y_1-b\right)-x_1 \cdot\left(k^{2}-1\right)}{k^{2}+1} x2=k2+12⋅k⋅(y1−b)−x1⋅(k2−1)
y 2 = 2 ⋅ k ⋅ ( b + x 1 k ) + y 1 ⋅ ( k 2 − 1 ) k 2 + 1 y_2=\frac{2 \cdot k \cdot\left(b+x_1 k\right)+y_1 \cdot\left(k^{2}-1\right)}{k^{2}+1} y2=k2+12⋅k⋅(b+x1k)+y1⋅(k2−1)
由于本问题中对称轴
y
=
x
+
1
y = x + 1
y=x+1特殊,实际计算对称点的代码很简单:
endc0=(end[1]-1)
endc1=(1+end[0])
Python代码
# coding=utf-8
import matplotlib.pyplot as plt
import numpy as np
from matplotlib.ticker import MultipleLocator
def plt_arrow(x_begin,y_begin,x_end,y_end,color,plt):#画带箭头的线段
plt.arrow(x_begin, y_begin, x_end - x_begin, y_end - y_begin,
length_includes_head=True, # 增加的长度包含箭头部分
head_width = 0.15, head_length =0.15, fc=color, ec=color)
def route1(path, plt,color):#处理一条路径,就是画一个子图
print(path)
lpath =len(path)
isCantelan = True
beginc0=0#保存越过限制线后箭头的起点坐标
beginc1=0
for i in range(lpath-1):
begin=path[i]
end=path[i+1]
#print(begin,end)
if(isCantelan):
if(begin[1]>begin[0]):#该点以后越过限制线:y=x,换颜色
isCantelan = False
plt_arrow(begin[0], begin[1], end[0], end[1], "orange",plt)
#求end点关于y=x+1的对称点
#print(end)
endc0=(end[1]-1)
endc1=(1+end[0])
#print("映射对称点坐标:",endc0,endc1)
plt_arrow(begin[0], begin[1], endc0, endc1, "magenta",plt)#画映射后的路径
#保存对称点坐标,为下一个箭头的起点
beginc0=endc0
beginc1=endc1
else:
plt_arrow(begin[0], begin[1], end[0], end[1], color,plt)#未越过限制线的路径
else:
plt_arrow(begin[0], begin[1], end[0], end[1], "orange",plt)
#求end点关于y=x+1的对称点
#print(end)
endc0=(end[1]-1)
endc1=(1+end[0])
#print("映射对称点坐标:",endc0,endc1)
plt_arrow(beginc0, beginc1, endc0, endc1, "magenta",plt)#画映射后的路径
#保存对称点坐标,为下一个箭头的起点
beginc0=endc0
beginc1=endc1
return isCantelan #返回该路径是否属于卡特兰路径
def buildpath(n):
#构造n为参数的全部路径
allpath=[]#全部路径
p =[]
p.append([0,0])
needstep=n
allpath.append(p)
for i in range(4):
newallpath=[]#新全部路径
for j,path in enumerate(allpath):
#print("\n当前路径:",j+1,path)
endpoint=path[-1][:]#取已有最后坐标作为下一步的基础
pxstep=endpoint[0]
pystep=endpoint[1]
if pxstep<needstep:
newpointx=endpoint[:]
newpointx[0]=newpointx[0]+1
pathx=path[:]
pathx.append(newpointx)
#print("扩展后路径:",pathx)
newallpath.append(pathx)
if pystep<needstep:
newpointy=endpoint[:]
newpointy[1]=newpointy[1]+1
pathy=path[:]
pathy.append(newpointy)
#print("扩展后路径:",pathy)
newallpath.append(pathy)
#print("\n第",i+1,"步所有路径",len(newallpath),"条:")
#[print(x) for x in newallpath]
allpath=newallpath[:]
return allpath
def drawallpath(allpath):
fig,axes = plt.subplots(2,3,figsize=(4.4,4),sharex=True,sharey=True)
plt.subplots_adjust(wspace=0.2,hspace=0.2)
xmajorLocator = MultipleLocator(1) #将x主刻度标签设置为1的倍数
ymajorLocator = MultipleLocator(1) #将y轴主刻度标签设置为1的倍数
for i in range(2):
for j in range(3):
id=i*3+j
#print(id)
axes[i,j].set_xlim(-0.1, 2.1)#指定画图范围
axes[i,j].set_ylim(-0.1, 3.1)
axes[i,j].xaxis.set_major_locator(xmajorLocator)
axes[i,j].yaxis.set_major_locator(ymajorLocator)
axes[i,j].xaxis.grid(True, which='major') #x坐标轴的网格使用主刻度
axes[i,j].yaxis.grid(True, which='major') #y坐标轴的网格使用次刻度
axes[i,j].grid(ls=":",c='gray',)#打开坐标网格
axes[i,j].spines['right'].set_visible(False)#去掉边框
axes[i,j].spines['top'].set_visible(False)
axes[i,j].spines['left'].set_visible(False)
axes[i,j].spines['bottom'].set_visible(False)
axes[i,j].scatter(2, 2, s = 30, color = 'r',marker='o')#目标终点
axes[i,j].scatter(1, 3, s = 30, color = 'purple',marker='o')#映射后的终点
lc=[0,3]#限制线坐标
axes[i,j].plot(lc,lc,color='blue', linewidth=0.7, linestyle='dashed') # 限制线
scx=[0,2]# 翻转对称线x坐标
scy=[1,3]# 翻转对称线y坐标
axes[i,j].plot(scx,scy,color='cyan', linewidth=0.7, linestyle='dashed') # 翻转对称线
isCantelan=route1(allpath[id],axes[i,j],"green")
if isCantelan:
axes[i,j].set_title(str(id+1)+' Cantelan')
else:
axes[i,j].set_title(str(id+1)+' Not')
plt.savefig(r"data\Catalan.png",dpi=72) #保存分辨率为72的图片
allpath=buildpath(2)
drawallpath(allpath)