はじめに
ほんの一時期、Quora というサイトをよくサーフィンしていたことがあって、その際に以下のような質問を見つけました。
集合Pの要素の重複あり組合せの各総和で表すことのできない最大の数を求めるプログラムをpythonでかきたいです。どなたか教えていただけないでしょうか？（引用元）
一応、この質問には回答（解答ではない）がついていますが、回答者が質問の意図を掴めず有耶無耶になっていました。確かに、質問の仕方というか言葉づかいが適切ではないところも多少あり、正確に回答できないと言う回答者の気持ちも分かりますが、恐らく次のような問題を Python で解きたかったんじゃないかと思います。
整数の集合 $P$ から重複ありで取り出した整数（要素）の総和として表すことのできない最大の整数を求めよ。
もう少し数学チックに表現すれば、次のように書き直せます。
整数の集合 $P=\{a_1,\ldots,a_n\}$ に対する線形和
S = S(a_1,\ldots,a_n) = \sum_{i=1}^n a_ix_i \quad
(x_i \in \mathbb{Z}_{\ge 0})
として表すことのできない最大の整数 $M=M(a_1,\ldots,a_n)$ を求めよ。
回答への質問者の返信を読む限り、最終的には質問者は $M(4,7,11)$ を求めるにはどのようなコードを書けばよいかが知りたいんだと思います。
このとき、質問者の返信の中にある
集合の一番小さい要素の回数（この場合は4）連続で数字をあらわすことができればその後の数字は表すことが可能
という一文も理解できます。つまり、連続する 4 個の整数が $S=4x_1+7x_2+11x_3$ の形で表せた場合、それらに 4 を加えた整数も $S$ の形で表せることに変わりはないので、それ以降のすべての整数が $S$ の形で表せることになります。
ちなみに $M(4,7,11)=17$ です。言い換えれば、18 以上の整数は 4, 7, 11 の線形和として必ず表せます。さらに言えば、$11=4+7$ なので 11 が有っても無くても結果は同じです。
2 変数の場合
表記をそれぞれ $(a_1,a_2,x_1,x_2)\to(a,b,x,y)$ に書き換えます。すなわち、$S=ax+by$ とします。
また、$a$ や $b$ が負の場合を考え出すと面倒だったので、ひとまずはどちらも 1 以上の整数である場合のみを考えます。
少なくとも一方が 1 の場合
$a$ と $b$ の少なくとも一方が 1 の場合、1 でないほう（両方とも 1 のときはいずれか一方）の変数を 0 に固定すれば、$S$ は 0 以上のすべての整数を表すことができるので $M=-1$ となります。
以降では、$a$ と $b$ はともに 2 以上であるとします。
最大公約数が 2 以上の場合
$a$ と $b$ の最大公約数を $d\ge 2$ とすると、$S$ は $d$ の倍数です。つまり、$d$ の倍数以外の整数は表すことができないので、$M$ は存在しません（$M=\infty$ になってしまう）。
上記以外の場合
$a$ と $b$ はともに 2 以上で、互いに素であるとします。このとき、Bézout の補題により $ax_0+by_0=1$ を満たす整数 $x_0,y_0$ が存在します。特に、任意の整数 $N,k$ に対して
N = a(Nx_0 - bk) + b(Ny_0 + ak)
\tag{1}
が成り立ちます。右辺は $a$ と $b$ の線形和になっているので、これの括弧内がともに 0 以上となる $k$ が存在しない $N$ のうち最大のものが $M$ です。
そのような $M$ が必ず存在すること、つまり、ある整数以上のすべての整数は $S$ の形で表せることを示します。
$a,b\ge 2$ および $ax_0+by_0=1$ より $x_0y_0\lt 0$（正負が異なる）です。また、$N=1$ のもとで $k$ を適切に選び、式（1）の括弧内を改めて $x_0,y_0$ と置き換えれば $x_0\gt 0$（つまり $y_0\lt 0$）にできます。
このとき、式（1）の括弧内がともに 0 以上となるためには、以下の不等式が成り立てば OK です。
- $Nx_0-bk\ge 0$
- $Ny_0+ak\ge 0$
これを $k$ について解くと $-\frac{Ny_0}{a}\le k\le\frac{Nx_0}{b}$ となります。これを満たす整数 $k$ が必ず存在するためには、両辺の差が 1 以上であればよいので、
\frac{Nx_0}{b} - \left(-\frac{Ny_0}{a}\right) = \frac{N}{ab} \ge 1
より $N\ge ab$ となります。すなわち、$ab$ 以上の任意の整数 $N$ に対しては式（1）の括弧内がともに 0 以上となる整数 $k$ が存在するので、$N$ は $S$ の形で表せます。
一方、少なくとも整数 1 は $S$ の形で表せません。
従って、$S$ の形で表せない最大の整数 $M$ は 1 以上 $ab$ 未満の間に必ず存在します。
ただし、これ以上 $M$ を特定または範囲を絞ることは難しそうです。実際、$a$ と $b$ の値によって $M=1$ のときもあれば $M=ab-1$ のときもあります。
M を求めるコード
import math
import sympy as sp
def find_Mab(a, b):
assert a >= 1 and b >= 1
# 少なくとも一方が 1
if a == 1 or b == 1:
return -1
# 最大公約数が 2 以上
if math.gcd(a, b) >= 2:
return None
# 上記以外
x0, y0, _ = sp.gcdex(a, b)
for N in range(a * b - 1, 1, -1):
k_lower = math.ceil(-N * y0 / a)
k_upper = math.floor(N * x0 / b)
if k_upper < k_lower:
return N
return 1
find_Mab(9, 7)
# 47
find_Mab(4, 7)
# 17
一般（多変数）の場合
恐らく一般的な解法はなさそうです。ここでは素直に以下の手順で求めていきます。なお、2 変数の場合と同様、各 $a_i$ は 1 以上であるとしておきます。
- 集合 $P$ の中に 1 が含まれている場合、$M=-1$ を出力する。
- 相異なる 2 つの整数 $a_i,a_j\in P$ に対して上記の 2 変数の場合を適用し、$M_{ij}$ を求める。
- 各 $M_{ij}$ の中で最小のものを $M_\min$ とし、$M_\min$ から 2 までの整数を順に $S$ の形で表せないか調べる。
- 表せないものが見つかればそれを $M$ として出力し、見つからなければ $M=1$ として出力する。
import itertools as it
def find_M(*A):
A = sorted(set(A))
for a in A:
assert a >= 1
# 少なくとも一つが 1
if 1 in A:
return -1
# 要素数が 1
if len(A) == 1:
return None
# 要素数が 2
if len(A) == 2:
return find_Mab(*A)
# それ以外
M_min = math.inf
for ai, aj in it.combinations(A, 2):
if (ai > M_min) or (aj > M_min):
continue
Mij = find_Mab(ai, aj)
if Mij is None:
continue
if Mij < M_min:
M_min = Mij
# 最大公約数が 2 以上
if M_min == math.inf:
return None
# M_min 以上の整数は不要なので除外
A = [ a for a in A if a <= M_min ]
# 組み合わせを計算
for N in range(M_min, 1, -1):
if not exist_combination(N, A):
return N
return 1
def exist_combination(N, A):
if len(A) == 1:
return N % A[0] == 0
A_sub = A[:-1]
an = A[-1]
for k in range(0, N + 1, an):
L = N - k
if L == 0:
return True
if exist_combination(L, A_sub):
return True
return False
find_M(4, 7, 11)
# 17
find_M(12, 35, 51, 70)
# 163
おわりに
もし間違いがあったり、より効率の良い方法があればぜひ教えてください。