Apache2.0 ライセンスとなった数理最適化ソルバー SCIP と Cbc を手持ちの問題で比較する
今年の11月に開催された数理最適化ソルバーの SCIP の ワークショップ において、 SCIPのライセンスが Apache 2.0 になることが発表されました。 SCIP は元々 ZIB という独自のライセンスで、学術利用でのみフリーで利用可能でしたが、今回のライセンス改定により、無償での商用利用も可能となります。アナウンスのあった11月の時点では、あくまで次期バージョンからの改訂という話だったのですが、先日その 8.0.3 がリリースされ、名実ともにApache 2.0 ライセンスとなったようです。
SCIP の他にもオープンソースな数理最適化ソルバーはいくつかあります[1]が、商用のもの[2]と比較すると、お世辞にも性能が良いとは言えない状況でした。そのため、解ける問題のサイズも限られたものになってしまい、数理最適化が普及しにくい要因の一つとなっていました。今回のライセンス改定により、オープンソースの範囲内で解ける問題サイズが増えそうです。
そこで、本記事では、最近自分が定式化して遊んだ問題を題材にして、どれくらい性能が上がるのかを見てみたいと思います。
SCIP のベンチマーク性能
独自の問題を解いて比較する前に、標準的なベンチマーク問題で複数のソルバーの性能を比較した結果を見てみましょう。
Cbc と比較して SCIP は頑張っていますが、Gurobi とは大きな差があることがわかります。今回目指せるのもこの辺りになりそうです。
SCIP のインストール
SCIP のバイナリは、ダウンロードページ から入手できます。
ただし、以下の記述がある通り、環境が Mac の場合は libgfortran などを用意する必要があります。
For some of the mac packages you need the libgfortran and libtbb libraries, and bison. You can install them for example via homebrew:
brew install bison
,brew install gcc
andbrew install tbb
before installing the scipoptsuite.
libgfortran は gcc パッケージに含まれていますが、対応している gcc のバージョンは 7 か 10 と少し古いものになっていて、M1 Mac で brew
を使ってインストールしようとすると、エラーとなってしまいます。
そこで、今回はソースコードからインストールすることにしました。といっても、インストール自体はとても簡単で、以下のコードで一発でした
# Boostがないとおこられる
brew install boost
# SOPLEXのインストール
git clone https://github.com/scipopt/soplex.git
cd soplex
mkdir build
cd build
cmake ..
make
sudo make install
# SCIPのソースコードをclone(ここでは master ブランチを利用)
cd ../..
git clone https://github.com/scipopt/scip.git
# ビルドディレクトリ作成
cd scip
mkdir build
cd build
# ビルド. DAUTOBUILD=ON で足りない独自ライブラリを探してダウンロードしてくれる
cmake .. -DAUTOBUILD=ON
make
# システムにインストール
sudo make install
インストール後に別ターミナルを立ち上げて scip
と打ち込み、プロンプトが表示されることが確認できれば、インストールは完了です。
比較1: 披露宴座席配置問題
数理最適化 アドベントカレンダー 2022 9日目の記事で、披露宴の座席配置を数理最適化を用いて決定した話が紹介されていました。その記事の中で高々20人であっても数理最適化ソルバーを使って厳密解を見つけるのは難しい点に言及されていました。
その話を受けて、数理最適化ソルバーを使ってどれくらいの問題を解くことができるか、現実的にはどういう戦法を取るのかという記事も公開されました。
実は自分も「自分だったらどうするかなぁ」ということを考えて定式化していたので、この問題を取り上げたいと思います。
解きたい問題
座席間の近さ
定式化
定式化は上記の記事のものとよく似ていますが[3]、論理積と
0 | 0 | 0 | 0 |
0 | 1 | 0 | 0 |
1 | 0 | 0 | 0 |
1 | 1 | 1 | 1 |
定数
-
:S_{i,j} とi の知り合い度合いj -
: 座席C_{a,b} とa の席の近さb
変数
-
x_{i,a} \in \{0, 1\} - ゲスト
が座席i にアサインされたときにa となるフラグ1
- ゲスト
-
z_{i,j,a,b} \in \{0, 1\},\ \mathrm{for}\ i < j - ゲスト
が座席i 、ゲストa が座席j にアサインされたときにのみb となる補助変数1 - 言い換えると
z_{i,j,a,b} = \min(x_{a,i},x_{j,b})
- ゲスト
目的関数
ゲスト
\mathrm{maximize} \sum_{i,j,a,b}z_{i,j,a,b}S_{i,j}C_{a,b}
制約条件
-
(各ゲストに1つ座席が与えられる)\sum_{a} x_{i,a} = 1 -
(全ての座席がうまる)\sum_{i} x_{i,a} = 1 -
(z_{i,j,a,b} \leq x_{i, a} がx_{i,a} のとき、0 もz_{i,j,a,b} )0 -
(z_{i,j,a,b} \leq x_{j, b} がx_{j,b} のとき、0 もz_{i,j,a,b} )0
上記に加えて、
実装
model = pulp.LpProblem('wedding_seats_assignment', sense=pulp.const.LpMaximize)
x = np.array(
pulp.LpVariable.matrix('x', (guest_names, seats), cat=pulp.LpBinary),
)
z = np.array(
pulp.LpVariable.matrix(
'z', (guest_names, guest_names, seats, seats), cat=pulp.LpBinary)
)
for i in range(x.shape[0]):
for j in range(x.shape[0]):
if i >= j:
z[i, j] = 0
z = np.tril(z.T).T
np.fill_diagonal(z, 0)
# 目的関数の設定
z_s = np.tensordot(guest_distances, z, axes=([0, 1], [0, 1]))
z_s_c = np.tensordot(z_s, seat_distances, axes=([0, 1], [0, 1])).sum()
model += z_s_c
# 制約条件の追加
for i in range(x.shape[0]):
model += pulp.lpSum(x[i, :]) == 1
for a in range(x.shape[1]):
model += pulp.lpSum(x[:, a]) == 1
for i in range(x.shape[0]):
for j in range(x.shape[0]):
if i >= j:
continue
for a in range(x.shape[1]):
for b in range(x.shape[1]):
model += z[i, j, a, b] <= x[i, a]
model += z[i, j, a, b] <= x[j, b]
model.solve(SOLVER)
比較結果
前述の通りゲストが20人いると単純に解くことはできないので、ゲスト数
Cbc | 157.53sec | 233.13sec | 711.03sec | > 6000sec |
SCIP | 10.23sec | 12.82sec | 27.22sec | 166.24sec |
ざっくりと5倍から10倍10倍以上高速化していそうです[3:1]。また、この手の問題はサイズが大きくなると急激に時間がかかるようになるものですが、その「急激に時間がかかるようになる」サイズも SCIP のほうが大きそうです。
なお Cbcでは
比較2: 最長しりとり問題
今年のアドベントカレンダーでさだまさしさんの曲タイトルでしりとりをしました。
上記の記事のポテンシャル最適化のコードをSCIPで動かすとどうなるかをみてみようと思います
解きたい問題
さださんの楽曲(ネット上でタイトルを確認できた 572曲)のタイトルでしりとりをした場合に最大でどれくらい続けられるかを求める
定式化
以下のように定式化しました。詳細は記事をご確認ください。
定数
-
: タイトル数N\in \mathbb{N}
変数
-
: タイトルx_{ij} \in \{0, 1\} からタイトルi につなげる時はj 、つなげない場合は1 となる変数0 -
: タイトルz_{i} \in \{0,\dots,N\} の訪問順i
目的関数
曲と曲がつながった数、つまり
\mathrm{maximize} \sum_{ij}{x_{ij}}
制約条件
-
(各タイトルは最大1回だけ利用可能)\sum_{J} x_{JI} \leq 1 \quad \mathrm{for}\ \forall I -
(あるタイトルにつづくタイトルは1つだけ)\sum_{J}x_{IJ} \leq 1 \quad \mathrm{for}\ \forall I -
(部分巡回路除去)z_{I} + 1 - N(1 - x_{Ij}) \leq z_{j}\ \forall I,\ j -
(しりとりにならない経路は利用できない)x_{i,j} = 0 \ \forall (i, j) \notin C
比較結果
こちらもタイトル数を変えて実行時間を測ってみました
Cbc | 0.38sec | 3.83sec | 54.34sec | 46.06sec |
SCIP | 0.21sec | 2.48sec | 142.59sec | 366.43sec |
面白いことに、
まとめ
Apache 2.0 ライセンスになった SCIP のインストール方法の紹介と Cbc を、過去自分が書いたコードベースでの比較結果をまとめました
- SCIP のインストールはとても簡単でお手軽。PuLP であれば、そのままソルバーを切り替えられるので、お手軽に高速化できます
- ベンチマークを見る限り、Cbc よりも SCIP の方が高速です。ただし、すべての問題において一律に高速というわけではなないこともわかりました。PuLP を使っている限りではソルバーの切り替えは簡単なので、問題毎に比較して良さそうなものを選択するとよさそうです
Discussion