python 组合贷 python组合数计算cmn=n!

文章目录

• 一、组合相关运算

• 1.1 组合数
• 1.2 伯努利数
• 1.3 卡特兰数（Catalan number）
• 1.4 欧拉数
• 1.5 第一类 Stirling 数
• 1.6 第二类 Stirling 数
• 1.7 Bell 数
• 1.8 整数的拆分

说明：本文延续上一篇文章，继续探索Python 3.X 中的数论计算模块 NZMATH 。

一、组合相关运算

1.1 组合数

组合数形如 $C_{n}^{m} = \frac{n!}{(n-m)!m!}$ ，计算 $C_{987}^{567}$

from nzmath import *
combinatorial.binomial(987,567)

输出结果为：$\begin{aligned} &56790077601071954516435076674439203693823757066221415014227080186952412759 \\ &74709993113491804414559512571055298475035368473091725522265486304592753587 \\ &64800433890809069149092493232964111334092617959312972975337288228379148827 \\ &014083315194945568254182774638098968556205453536200706386074651328512 \end{aligned}$ 关于 $100!$

combinatorial.factorial(100)

我们能进行大阶乘运算，比如还可以计算 $1000!$、$10000!$、$100000!$

CombinationIndexGenerator 为组合数的迭代器，组合数的情况罗列如下

list(combinatorial.combinationIndexGenerator(5,3))

输出结果为

[[0, 1, 2],
 [0, 1, 3],
 [0, 1, 4],
 [0, 2, 3],
 [0, 2, 4],
 [0, 3, 4],
 [1, 2, 3],
 [1, 2, 4],
 [1, 3, 4],
 [2, 3, 4]]

PermutationGenerator 为全排列数的迭代器，全排列数 $A^3_3$

list(combinatorial.permutationGenerator(3))

输出结果为

[[0, 1, 2], [0, 2, 1], [1, 0, 2], [1, 2, 0], [2, 0, 1], [2, 1, 0]]

1.2 伯努利数

伯努利数可以由微分式定义 $B_n = \frac{d^n}{dx^n}\left.\left( \frac{x}{e^x-1} \right) \right|_{x=0}$ 其中，$B_0=1$ ，$B_1=-\frac{1}{2}$ ，$B_2=\frac{1}{6}$ ，$B_3=0$ ，$B_4=-\frac{1}{30}$ 、$B_5=0$ 。值得一提的是伯努利数与 Zeta 函数的关系 $\zeta(2n) = \frac{(-1)^{n-1}(2\pi)^{2n}}{2(2n)!} B_{2n}$ 另一个总所周知的结论是，对于 $n\geq1$ ，有 $B_{2n+1}=0$

现在计算第 $100$

combinatorial.bernoulli(100)

输出结果为：

Rational(-36186471006482321128123819673322305121771572462967870388929762640041513818242476739343599169377521437613732462774398901923, 12748833821603357261103218279966734353954750)

其中，Rational(m, n) 表示分数 $\frac{m}{n}$

1.3 卡特兰数（Catalan number）

卡特兰数定义如下 $C_n = \frac{C_{2n}^n}{1+n} ,\quad n\geq 0$

combinatorial.catalan(200)

结果如下：

5122014932110168435825813115249447769376062806324872654056736267065921305612196023
99612226723109479418508961937154088

1.4 欧拉数

根据 $\mathrm{sech\,}(x)$ 在 $x=0$ 的泰勒展开即可引入欧拉数 $E_n$：$\mathrm{sech\,}(x) = \sum^{\infty}_{n=0} \frac{(-1)^n}{(2n)!} E_n x^{2n}$

combinatorial.euler(200)

输出结果为：

5995471635010873724785738304825045455122720583803852868566645394518582604129494754
8961407177543691976315500802573962474205261301751817830467481717661644330627229689
8776839142684308452800807505458856711470377113494564236618553260087345930603179179
0252681485922762322793654848981287557319739693929515377915693818019564833879663352
33673223

1.5 第一类 Stirling 数

第一类斯特林数表示把 $n$ 个不同元素构成 $m$ 个圆的排列方案数，也是无符号的第一类斯特林数，记作 $S_u(n,m)$ 。带符号的第一类斯特林数可以简记为 $S_s(n,m)=(-1)^{n+m}S_u(n,m)$

for i in range(11):
    print(combinatorial.stirling1(10, i), end=' ')

输出结果为：

0 -362880 1026576 -1172700 723680 -269325 63273 -9450 870 -45 1

1.6 第二类 Stirling 数

第二类 Stirling 数涉及集合的拆分，表示将 $n$ 个不同的元素拆分成 $m$ 个集合的方案数，记为 $S(n,m)$ 。与第一类斯特林数不同的是，集合里面的元素是无序的，而圆排列是有序的。我们来计算 $S(100,49)$

combinatorial.stirling2(100,49)

得到的结果是：

525178013011792678224299033717178581734767648826742619698978792854110316239515066
9196100844267828152500

循环 $n=10$

for i in range(11):
    print(combinatorial.stirling2(10,i), end=' ')

输出结果为：

0 1 511 9330 34105 42525 22827 5880 750 45 1

1.7 Bell 数

贝尔数的定义为 $b(n) = \sum^{n}_{i=0} S(n,i)$ 其中，$S(n,i)$

combinatorial.bell(200)

算得第 200 个贝尔数为

624748477619370179475169176526172040868676631541156088557319107233592458655068693
596197197036539240046478541880130367251305801992731793116433134463929399319082714
196961652168013602803041772182448739785880587721883539643696750170168607678730423
675712080858884633763266376601388

1.8 整数的拆分

将一个整数 $n$ 拆分成不同的和式，我们把拆分的和式数记为 $p(n)$ ，比如 $3=1+2=1+1+1$ ，故 $p(3)=3$ ，显然这是一个关于 $n$ 的无序拆分。现求 $n=20000$ 的拆分数 $p(20000)$

combinatorial.partition_number(20000)

输出结果为：

252114813812529697916619533230470452281328949601811593436850314108034284423801564
956623970731689824369192324789351994903016411826230578166735959242113097

关于整数拆分，印度天才数学家拉马努金做出了卓越的成果：$p(n) \sim \frac{1}{4n\cdot \sqrt{3}} \exp\left(\pi \sqrt{\frac{2n}{3}}\right) , \; n\to \infty \tag{1}$ 根据上述公式，记 $\tilde{p}(n)$ 为右边的估计式，则我们可以计算误差百分比 $\frac{\tilde{p}(n)-p(n)}{p(n)}$

# 计算误差
def p(n):
    return (1/(4*n*np.sqrt(3)))*np.exp(np.pi*np.sqrt((2*n)/3))

rest = []
for i in range(10, 2001):
    rest.append((p(i)-combinatorial.partition_number(i))/combinatorial.partition_number(i))

# 绘图
import matplotlib.pyplot as plt

plt.figure(dpi=200)

x1 = [i for i in range(10, 2001)]
y1 = rest
plt.plot(x1, y1, color='b', linewidth=1)

x2 = [i for i in range(10, 2001, 50)]
y2 = [rest[i] for i in range(0, len(rest), 50)]
plt.scatter(x2, y2, color='m', s=10)

plt.show()

得到一条单调递减的相对误差曲线：

发现 $(1)$ 式的估计是很不错的。当 $n=10000$ 时，此时相对误差仅为 $0.004446$

最后，我们计算整数拆分的所有情况：

for part in combinatorial.partitionGenerator(6):
    print(part)

输出结果为：

(6,)
(5, 1)
(4, 2)
(4, 1, 1)
(3, 3)
(3, 2, 1)
(3, 1, 1, 1)
(2, 2, 2)
(2, 2, 1, 1)
(2, 1, 1, 1, 1)
(1, 1, 1, 1, 1, 1)

文章暂时先写到此处吧，欢迎读者浏览本系列文章的后续内容。