③整数線形最適化問題 トロピカル幾何学の勉強
まだ生成AIでトロピカル演算が得意なのが無いが、意味は分かった。
なるほどね、確かに幾何的に解くんだな。イメージできた。
ひょっとするとClaudeは直接に解いてくれちゃうかもしれないが。
どうするか、ゲーム理論の時みたいにMathematicaのコードを書いてもらおうかな。
P10の例1.2.4
トロピカル演算で整数線形計画法を解きます。
Au=bの元でwuを最小化したい。
A=(a_ij),b=(b_1,...,b_d)^Tを非負整数を成分に持つ列ベクトルとする。
目的は非負整数列ベクトルu=(u_1,...,u_n)^Tで、
次の最適問題の解となるものである。
目的：Au=bの条件のもとでwuを最小にする。
行列Aのすべての列の成分の和が同じ数αに等しく、
b_1+...+b_d=mα u_1+...+u_n=mと仮定する。
意味はこうだ。
頑張ったが、どの生成AIでもトロピカル幾何学を用いて直接には解いてくれなかった。でも、解説のおかげで問題のイメージは掴めた。
トロピカル演算を使った整数線形計画問題の解法について考えていきましょう。まず問題を整理します。
与えられた条件:
A = (a_ij)は非負整数行列
b = (b_1,...,b_d)^T は非負整数ベクトル
行列Aの各列の成分の和 = α (定数)
Σb_i = mα
Σu_i = m
求めるu = (u_1,...,u_n)^T は非負整数ベクトル
目的:
Au = b を満たし、w*u を最小化する
このような問題に対して、
トロピカル演算のアプローチを使用すると以下の様に解くことができます:
トロピカル半環での演算を定義:
⊕ (トロピカル加算): min演算
⊗ (トロピカル乗算): 通常の加算
目的関数の変換: w*u をトロピカル演算で表現すると:
w_1⊗u_1 ⊕ w_2⊗u_2 ⊕ ... ⊕ w_n⊗u_n
制約条件の変換:
Au = b を各成分ごとに:
min{a_i1 + u_1, a_i2 + u_2, ..., a_in + u_n} = b_i
解法のステップ:
a) 各変数u_iに対して可能な値の範囲を設定
b) 制約条件Σu_i = mを考慮
c) トロピカル演算の性質を利用して、最適解を探索
この問題の特徴として、行列Aの各列の和が等しいことと、
Σb_i = mα
という条件があります。これらの条件は、解空間を大きく制限し、効率的な解の探索を可能にします。
これは意図に反して初等的に解いてるけど意味は分かりやすい。
では、具体的な数値で解いていきましょう。
まず、与えられた値を確認します：
A = [4 3 2 1 0] [0 1 2 3 4]
b = [5] [7]
w = (2,5,11,7,3)
問題の条件を確認：
各列の和 α = 4 (全ての列で4)
b_1 + b_2 = 12 = 3α なので、m = 3
したがって u_1 + u_2 + u_3 + u_4 + u_5 = 3
制約方程式：
4u_1 + 3u_2 + 2u_3 + u_4 + 0u_5 = 5 0u_1 + u_2 + 2u_3 + 3u_4 + 4u_5 = 7
非負整数解を探索するために、以下のステップで解きます：
u_1 + u_2 + u_3 + u_4 + u_5 = 3
という条件と、
各u_iが非負整数という条件から、可能な組み合わせを列挙します。
これらの組み合わせの中から、制約方程式を満たすものを見つけます。
その中から w*u が最小となる解を選びます。
可能な解を探索すると：
u = (0,0,1,1,1) のとき：
4*0 + 3*0 + 2*1 + 1*1 + 0*1 = 3 = 5 (×)
0*0 + 1*0 + 2*1 + 3*1 + 4*1 = 9 ≠ 7 (×)
u = (1,0,0,1,1) のとき：
4*1 + 3*0 + 2*0 + 1*1 + 0*1 = 5 (○)
0*1 + 1*0 + 2*0 + 3*1 + 4*1 = 7 (○)
これは制約を満たします。他の組み合わせも確認すると、これが有効な解の一つとなります。
このとき、wu = 2*1 + 5*0 + 11*0 + 7*1 + 3*1 = 12
他の可能な組み合わせも確認しましたが、この解が最適解となります。
したがって、最適解は： u* = (1,0,0,1,1) 最小値：wu = 12
この解は：Au = b を満たす
各成分が非負整数
成分の和が3
w*u を最小化する
という全ての条件を満たしています。
これはつまり、今まではPythonの特殊なパッケージ（PuLPとか）を入れないと解けなかった問題が、そのままdefで定義したルーチンだけで動く事になる。実際、現場のPCだとpipが通じなくてかなり面倒な事があって…
Mathematicaで解くのも少し癖が有る。
トロピカル演算です。 p=2x_1^4⊕5x_1^3x_2⊕11x_1^2x_2^2⊕7x_1x_2^3⊕3x_2^4です。 p⊗p⊗pを求めるMathematicaのコード
同じ次数の項目は、通常演算で足して表示。
多項式を文字列ではなく、Mathematicaネイティブのデータで処理して。
12x1^5x2^7 で x1^5x2^7 の係数が 12 なので
答えは12
(*トロピカル演算の定義*)tropicalPlus[x_, y_] := Min[x, y]
tropicalTimes[x_, y_] := x + y
(*初期多項式の定義*)
p = 2 x[1]^4 + 5 x[1]^3 x[2] + 11 x[1]^2 x[2]^2 + 7 x[1] x[2]^3 +
3 x[2]^4
(*多項式からリスト形式への変換関数*)
toTermList[poly_] := CoefficientRules[poly, {x[1], x[2]}]
(*リスト形式から多項式への変換関数*)
toPoly[rules_] :=
Sum[rules[[i, 2]] x[1]^(rules[[i, 1, 1]]) x[
2]^(rules[[i, 1, 2]]), {i, Length[rules]}]
(*トロピカル乗算の実装*)
tropicalProduct[p1_, p2_] :=
Module[{terms1 = toTermList[p1], terms2 = toTermList[p2],
result},(*全ての項の組み合わせに対してトロピカル乗算を適用*)
result =
Flatten[Table[{terms1[[i, 1]] + terms2[[j, 1]],
tropicalTimes[terms1[[i, 2]], terms2[[j, 2]]]}, {i,
Length[terms1]}, {j, Length[terms2]}], 1];
(*同じ指数を持つ項をグループ化して最小値を取る*)
result = GroupBy[result, First -> Last, Min];
(*結果を多項式形式に変換*)toPoly[Normal[result]]]
(*p\[CircleTimes]p の計算*)
p2 = tropicalProduct[p, p]
(*p\[CircleTimes]p\[CircleTimes]p の計算*)
p3 = tropicalProduct[p, p2]
(*結果の表示*)
Print["p\[CircleTimes]p\[CircleTimes]p の結果:"];
Print[p3]p\[CircleTimes]p\[CircleTimes]p の結果:
6 x[1]^12+9 x[1]^11 x[2]+12 x[1]^10 x[2]^2+11 x[1]^9 x[2]^3+7 x[1]^8 x[2]^4+10 x[1]^7 x[2]^5+13 x[1]^6 x[2]^6+12 x[1]^5 x[2]^7+8 x[1]^4 x[2]^8+11 x[1]^3 x[2]^9+17 x[1]^2 x[2]^10+13 x[1] x[2]^11+9 x[2]^12