やさしい母関数(Generating Function)(その5)で母関数をNumpyのconvolveを使って解きました。今回はその応用です。
色々な分割数
今回は色々な分割数の問題を母関数を用いて解いて行きたいと思います。
nをいくつかの 自然数の和で表す |
例 | |
---|---|---|
(1) | 制限なし OEIS-A000041 |
n=4, 5通り 4, 3+1, 2+2, 2+1+1, 1+1+1+1 |
(2) | 最大m個の和 |
n=4,m=2, 3通り 4 3+1, 2+2 |
(3) | 自然数の最大値がk |
n=4,k=2, 3通り 2+2, 2+1+1, 1+1+1+1 |
(4) | 奇数の和 OEIS-A000009 |
n=4, 2通り 3+1, 1+1+1+1 |
(5) | 互いに異なる自然数 OEIS-A000009 |
n=4, 2通り 4, 3+1 |
(6) | 素数の和 |
n=4, 1通り 2+2 |
(7) | 自己共軛 (self-conjugate) OEIS-A000700 |
n=8, 2通り![]() |
母関数(Generating Function)からPythonのコード
やさしい母関数(Generating Function)(その4)で示したように分割数の母関数は
\begin{align}
&\large \prod_{k=1}^{\infty}(\frac{1}{1-x^k}) \\
&= (1+x+x^2+\cdots) (1+x^2+x^4+\cdots) (1+x^3+x^6+\cdots)\cdots
\end{align}
やさしい母関数(Generating Function)(その5)で以下の1の「制限なし」のPythonのコードを紹介しました。
N = 4
p = [1]
for n in range(1,N+1):
p = np.convolve(p,[1 if i % n == 0 else 0 for i in range(N+1)])[:N+1]
print(f"Answer: {p[N]}")
# Answer: 5
色々な分割数を求めるコード
このコードをベースにして(1)(3)(4)(5)(6)を求めるように変更します。
- 分割するときに使える数をparts_inで指定できるようにする
- 互いに異なる数を使うフラゲdiscinctがTrueのときの処理を加える
import numpy as np
import sympy as sp
def partitions(N, parts_in=None, distinct=False):
if parts_in == None: parts_in = range(1,N+1)
p = [1]
for n in parts_in:
p1 = np.zeros(N+1, dtype="int64")
if distinct: p1[:n+1:n] = 1
else: p1[::n] = 1
p = np.convolve(p, p1)[:N+1]
return p
N, k = 4, 2
print(f"# N = {N}")
print(f"# 1: {partitions(N)[N]}通り")
print(f"# 3: {partitions(N, parts_in=range(1,k+1))[N]}通り")
print(f"# 4: {partitions(N, parts_in=range(1,N+1,2))[N]}通り")
print(f"# 5: {partitions(N, distinct=True)[N]}通り")
print(f"# 6: {partitions(N, parts_in=list(sp.primerange(2,N+1)))[N]}通り")
# N = 4
# 1: 5通り
# 3: 3通り
# 4: 2通り
# 5: 2通り
# 6: 1通り
いくつかの分割恒等式
1.「個数」と「最大数」
「個数」と「最大数」に関する性質 (高校数学の美しい物語)によると
「k個の正の整数に分割する方法の数」と「最大の数がkであるようないくつかの正の整数に分割する方法の数」は等しい。
従って(2)と(3)は$m=k$のとき等しくなるので、(2)は(3)のコードで求めることができます。
2. オイラーの分割恒等式
オイラーの分割恒等式によると、
「互いに異なる自然数に分割する方法の個数」と「奇数の自然数に分割する方法の個数」が等しい
なので$(4)=(5)$が成り立ちます
3.自己共軛と相異なる奇数
自然数の分割#フェラーズ図形によると、
自己共軛(self-conjugate)な分割の総数は相異なる奇数への分割の総数に等しい
なので、$(7)=(4) and (5)$で以下のコードで自己共軛(self-conjugate)分割の数が求まります
N = 8
print(f"# 7: {partitions(N, parts_in=range(1,N+1,2), distinct=True)[N]}通り")
# 7: 2通り
いくつかの例題
このコードを使ってProject Euerの問題を解いてみます。
【Project Euler】Problem 76: 整数の分割
N = 100
print(f"# Problem 76: {partitions(N)[N]-1}")
# Problem 76: 190569291
【Project Euler】Problem 77: 整数を素数に分割
N = 2
while partitions(N, parts_in=list(sp.primerange(2,N+1)))[N] < 5000:
N += 1
print(f"# Problem 77: {N}")
# Problem 77: 71
【Project Euler】Problem 31: コインの合計
N = 200
print(f"# Problem 31: {partitions(N, parts_in=[1,2,5,10,20,50,100,200] )[N]}")
# Problem 31: 73682
(開発環境:Google Colab)