数理最適化の応用例と 実験的解析 東京海洋大学 久保 幹雄
実験対象 施設配置問題(k-median, k-center) グラフ関連(グラフ分割,グラフ彩色) 巡回路(巡回セールスマン) ロットサイズ決定 非線形施設配置 スケジューリング
k-median問題 min-sum型の施設配置問題 顧客数 n=200,施設数 k=20 (顧客上から選択) Euclid距離,x,y座標は一様ランダム n=200, k=5 の例
定式化 I: 顧客の集合 J: 施設候補地点の集合 Big M 弱い定式化
弱い定式化での結果 n=200,k=20 Optimize a model with 401 Rows, 40200 Columns and 80400 NonZeros 中略 Explored 1445 nodes (63581 simplex iterations) in 67.08 seconds Thread count was 2 (of 2 available processors) Optimal solution found (tolerance 1.00e-04) Best objective 1.0180195861e+01, best bound 1.0179189780e+01, gap 0.0099% Opt.value= 10.1801958607
上界と下界の変化 (弱い定式化)
強い定式化での結果 Optimize a model with 40401 Rows, 40200 Columns and 160400 NonZeros 中略 Explored 0 nodes (1697 simplex iterations) in 3.33 seconds (分枝しないで終了!) Thread count was 2 (of 2 available processors) Optimal solution found (tolerance 1.00e-04) Best objective 1.0180195861e+01, best bound 1.0180195861e+01, gap 0.0% Opt.value= 10.1801958607
知見 Big Mを用いない強い定式化が望ましい. この程度の式なら,必要な式のみを切除平面として追加するような小細工は必要なし.
k-center問題 min-max型の施設配置問題 100顧客,10施設(顧客上から選択) Euclid距離,x,y座標は一様ランダム k-center (n=30,k=3) k-median (n=30,k=3)
定式化
上界と下界の変化(n=100,k=10) min-max型の目的関数には注意!
k-被覆問題 被覆されていない(距離がθより大きい)顧客数 =1 顧客 i が被覆されていない パラメータ:=1 顧客iから施設jへの距離が θ以下
k-被覆問題を用いた二分探索法 UB:=距離の最大値, LB:=0 while UB – LB >ε: θ= (UB+LB)/2 if k-被覆問題の最適値= 0 then UB = θ else LB = θ
実験結果 k= ceil( n/10 ) k-center min-max k-median k-center 2分探索
知見 k-center問題(min-max型の最適化問題)に対しては100点あたりで組合せ爆発が起きる. k-center問題に対しては,k-cover問題を用いた二分探索が有効
グラフ分割問題 無向グラフ G=(V,E)が与えられたとき,|L|=|R|を満たすVの分割(L,R)で,LとRの間の枝の本数を最小にするものを求める問題. 応用 VLSI 設計, プログラム分割 JohnsonらがSimulated Annealing法の実験のために採用した第1の問題. 16
グラフ分割問題 半分に分けたい! 仲良し! 17
グラフ分割問題 目的関数値=2 18
二次最適化定式化 凸二次でないので通常は解けない(がGurobiだと解ける) 0-1 変数を二乗して凸関数に変換する方が良い
Google 3D Graph x*(1-y)+(1-x)*y x^2+y^2-2*x*y
線形化定式化 変数 y を追加して線形化
二次錘定式化 二次最適化定式化 線形化定式化
知見 Gurobiは,非凸二次最適化でも0-1変数なら自動的に凸最適化に変換 線形化した方が高速 二次錘最適化も有効であるが,線形化が簡単で有効 100点を超えるような問題例には,適用すべきでない.
最大クリーク問題 (最大安定集合問題) 応用 1993年度DIMACSチャレンジ問題の一つ 無向グラフ G=(V,E)が与えられたとき,最大位数の安定集合 S (Vの部分集合で間に枝のないもの)を求める問題. 応用 プリント基板テスト,パターンマッチング,符号理論,考古学および生物学のデータ解析,グラフ彩色問題の近似アルゴリズムおよび下界,多面体論(Kellerの予想),故障診断,etc. 1993年度DIMACSチャレンジ問題の一つ 24
最大クリーク問題 なるべく大人数でピクニックに行きたい! 仲良し! 25
最大クリーク問題 26
最大安定集合問題 仲が悪い! 27
最大安定集合問題の定式化
グラフ彩色問題 無向グラフG=(V,E)が与えられたとき,Vの安定集合への分割で分割数最小のものを求める問題. 応用 時間割作成,スケジューリング,周波数割当,レジスタ配分,プリント基板テスト,パターンマッチグ. 1993年度DIMACSチャレンジ問題の一つ Johnsonらの実験的解析の第2の問題 30
グラフ彩色問題 ・・・ クラス分けしたい! 仲が悪い! 31
グラフ彩色問題 32
定式化 弱い定式化
グラフ彩色問題に対する実験 目的:解の対称性を調べる 点数 n=40,彩色数上限 Kmax=10 ランダムグラフ G(n,p=0.5) 彩色数 8 が最適値
上界と下界の変化(原定式化) 点数 n=40,彩色数上限 Kmax=10 Optimize a model with 3820 Rows, 410 Columns and 11740 NonZeros Explored 17149 nodes (3425130 simplex iterations) in 1321.63 seconds
定式化の改良 model.addSOS(1,変数リスト) 対称性の除去(番号の小さい方の色を優先して使う!) 特殊順序集合(SOS: Special Ordered Set) Type 1 (いずれか 1つの変数が正になることの宣言)の追加 model.addSOS(1,変数リスト)
上界と下界の変化(対称性除去) Optimize a model with 3829 Rows, 410 Columns and 11758 NonZeros Explored 4399 nodes (1013290 simplex iterations) in 384.53 seconds MIPFocus=2(最適性保証優先) で67秒,MIPFocus=3(下界優先) で70秒
上界と下界の変化(+SOS) Optimize a model with 3829 Rows, 410 Columns and 11758 NonZerosExplored 109 nodes (58792 simplex iterations) in 22.02 seconds MIPFocus=2(最適性保証優先) で65秒,MIPFocus=3(下界優先) で126秒
彩色数 K を固定した定式化 悪い枝の数を最小化 (0だとK彩色) 枝の両端点 i,j が同色だと 0-1変数 zij が 1
K固定問題を用いた二分探索 UB, LB := 彩色数 K の上界と下界 while UB – LB >1: K= [ (UB+LB)/2 ] if Kを固定した問題の最適値 = 0 then UB := K else LB := K
対称性除去 +SOS 原定式化 二分探索 対称性除去
知見 解の対称性のある問題は,下界の改善がしにくいので,分枝限定法では解きにくい 対称性を除く工夫を入れると多少は改善 SOS(Special Ordered Set)制約でさらに改善 彩色数に対する二分探索は有効
巡回セールスマン問題(TSP) すべての点をちょうど1回通る最短巡回路 切除平面法で85,900点の実際問題(対称TSP)の最適解
Miller-Tucker-Zemlinの定式化
上界と下界の変化 (80点,Euclid TSP) 強化した式 でないと... 1日まわして Out of Memory!
結果 Optimize a model with 6480 Rows, 6400 Columns and 37762 NonZeros 中略 Cutting planes: Gomory: 62 Implied bound: 470 MIR: 299 Zero half: 34 Explored 125799 nodes (2799697 simplex iterations) in 359.01 seconds Optimal solution found (tolerance 1.00e-04) Best objective 7.4532855108e+00, best bound 7.4525704995e+00, gap 0.0096% Opt.value= 7.45328551084
非対称巡回セールスマンのベンチマークに対する結果(時間) 多品種フロー 単一品種フロー 強化MTZ MTZ
非対称巡回セールスマンのベンチマークに対する結果(成功数) 多品種フロー MTZ 強化MTZ 単一品種フロー
知見 式を持ち上げ操作で強化すると,高速化され,大きな問題例が解けるようになる. 多品種フロー定式化は下界は強いが求解時間は悪い 単品種フロー定式化は安定して良い. 問題例によっては1000点まで解ける!
対称巡回セールスマンのベンチマークに対する結果(時間) 分枝カット法 II (分枝時にカット追加) 分枝カット法 I (整数解のとき カット追加) 分枝限定+ 切除平面法
対称巡回セールスマンのベンチマークに対する結果(成功数) 分枝カット法 II (分枝時にカット追加) 分枝カット法 I (整数解のとき カット追加) 分枝限定+ 切除平面法
知見 分枝時に切除平面を追加する分枝カット法は遅い. 整数解を発見したタイミングで切除平面を加える分枝カット法が推奨される. 分枝限定法で整数解を出した後で,切除平面を加える簡便法も悪くない.
多品目ロットサイズ決定問題 段取り費用と在庫費用のトレードオフを最適化する多期間生産計画 多品目で共通の資源を使う容量制約付き問題は,MIPソルバーには難問(と言われてきた) T=30期,P=24品目:Trigeiro, Thomas, McClain (1989年)の最大のベンチマーク
ロットサイズ決定問題 標準定式化のフローモデル 生産量(t) 在庫量(t) 在庫量(t-1) 期 t 需要量(t) 弱い定式化の原因 生産量(t)≦大きな数 “Large M” ×段取りの有無(t) 在庫量(t-1)+生産量(t)=需要量(t)+在庫量(t) 0-1変数
ロットサイズ決定問題 施設配置定式化のフローモデル s期に生産してt期まで在庫される量 期 s 期 t 需要量(t) s期に生産してt期まで在庫される量 ≦需要量(t)×段取りの有無(s) s期に生産してt期まで在庫される量 = 需要量(t)
アルゴリズムと定式化の関係 (特殊形に対する) 強い定式化 or 多項式時間アルゴリズム 多項式時間で破っている 妥当不等式の同定 追加制約+つなぎ制約 容量制約なしロットサイズ決定 Wagner-Whitin 動的計画アルゴリズム 施設配置定式化 最短路定式化 (S,l)不等式 1機械重み付き完了 時刻和最小化スケジューリング どん欲解法(Simth’s rule) Super-modular不等式
定式化のサイズと強さの比較 標準定式化 施設配置定式化 変数の数 変数の数 弱い定式化 強い定式化 =線形計画緩和が 整数多面体と一致 制約の数 制約の数 (S,l)不等式 切除平面 追加した 制約の数 n: 期数 強い定式化
上界と下界の変化(標準定式化) 1800秒で最適解;これ以上大きな問題例は無理!
上界と下界の変化(施設配置定式化) 40秒で最適解;T=100でも大丈夫!
切除平面 標準 施設配置
切除平面 標準 施設配置
知見 従来では難問と言われてきたロットサイズ決定問題でもある程度までは大丈夫 施設配置定式化を推奨 分枝カット法は悪い 実務的には「打ち切り分枝限定法」を推奨(緩和固定法と比べて簡単で高性能)
非線形施設配置問題 (色々な定式化の比較) 多重選択 非集約型凸結合 対数個の変数を用いた非集約型凸結合 集約型凸結合 対数個の変数を用いた集約型凸結合 SOS : タイプ2の特殊順序集合
集約型凸結合 対数個非集約型凸結合 対数集約型凸結合 非集約型凸結合 多重選択 sos
知見 単純なタイプ2の特殊順序集合 (SOS Type 2)を用いたものが最も良い. 次いで,対数個の変数を用いた集約型凸結合定式化,対数個の変数を用いた非集約型凸結合定式化,集約型凸結合定式化. 多重選択を使うケースをよくみるが,使うべきではない.
スケジューリング問題 1機械リリース時刻付き総納期遅れ最小化問題 比較する定式化 線形順序付け定式化 時刻添え字定式化 離接定式化 (Big M使用)
時刻添え字 線形順序付け 離接定式化
知見 比較的簡単な1機械問題でも200ジョブ程度で組合せ爆発 他のスケジューリング問題で実験したが,ベンチマーク問題例は全滅 現状では,スケジューリングは混合整数最適化では無理
Multi-objective Optimization Generate a set of non-inferior (Pareto) solutions by using Gurobi (that can find multiple solutions by one run of B&B) Make one objective by a linear combination of multiple objective functions Change objective functions to constraints and sliding the right-hand sides of the constraints Minimize the distance from ideal points
Non-inferior (Pareto) solutions Second Obj. Function Solutions in the region are dominated by solution A Linear Obj. Func. A Missed non-inferior point First Obj. Function
Best solution Under the time const. Generated solutions by the constraint sliding method (2-obj. TSP)
Constraint sliding (segmentation) and ideal point methods