1
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

ラグランジュ緩和で上界を作ってIntroduction to Heuristics Contestを解く

1
Posted at

ラグランジュ緩和で上界を作ってIntroduction to Heuristics Contestを解く

はじめに

以前、Introduction to Heuristics Contest(以下Intro)を題材に、AHC延長戦の伝説"saitodevel01"の機械学習手法を再現して1位を奪い返した話という記事を書きましたが、またしてもocchanというユーザにより、延長戦1位が更新されました。

2位.png

しかも、驚くべきことに、occhan氏の提出にはニューラルネットワーク(以降NNと表記)は入っていません。

この提出の軸は、「一部の制約をいったん外し、代わりに違反へ仮の価格を付ける方法」と、「最終スコアはこれ以下」という値を使うビームサーチです。前者を後ほどラグランジュ緩和、後者の値を上界という名前で説明します。

さらに、筆者によるその後の実験で、価格更新部分に安定化処理を追加し、少しだけスコアを伸ばすことができました。

1位.png

方法 延長戦スコア
筆者のラグランジュ緩和版 128,582,766
occhan氏の提出 128,561,553
筆者のNNビーム+後処理・探索効率化版 128,353,276
前回記事の最終提出 128,286,311

また、汎用テクニックとしてまとめられそうな部分をライブラリ化し、小問題を解く例からIntroで応用する流れを説明します。
ライブラリとその活用コードはサンプルコードとして公開しています。使えそうな問題があればご活用ください。

まずは小さな多次元0/1ナップサック

問題

品物ごとに価値、重さ、体積があります。それぞれの品物を選ぶか選ばないか決め、重量上限と体積上限を守りながら価値の合計を最大化します。

次の小さな問題を考えます。重量上限と体積上限は、どちらも10です。

品物 価値 重さ 体積
1 14 10 1
2 14 1 10
3 10 5 5
4 10 5 5

品物1は重量をほぼ全部使い、品物2は体積をほぼ全部使います。品物3と4は単体の価値では負けますが、両方選ぶと重量10、体積10、価値20になります。

ここで、単純な貪欲法を考えてみましょう。各品物について「重量容量の何割を使うか」と「体積容量の何割を使うか」を足し、価値をその値で割ります。その比が大きい品物から順に選ぶと、品物1を選んで14点になります。この方法では、品物3と4を選ぶ20点の解に届きません。そこで、上界を使ってよりよい選び方を探します。

上界と下界とは

この例の結果を見る前に、上界と下界という言葉を説明します。

この例では、後で説明するビームサーチが、制約を守る20点の選び方を実際に作ります。したがって、まだ本当の最適値が分からなくても、「最適値は少なくとも20点」と断言できます。このように、最適値がその値以上だと保証する値を、最大化問題の下界と呼びます。

一方、後で説明する計算により、「どのように品物を選んでも25.492点を超えない」と保証できます。このように、最適値がその値以下だと保証する値を上界と呼びます。

実際に作れた20点  <=  本当の最適値  <=  超えないと保証できる25.492点
       下界                                      上界

上界は、実際に作れる解のスコアとは限りません。「正解はここより上にはない」という目印です。下界は反対に、「この点数の解は実際に作れたので、正解は少なくともここまではある」という目印です。本記事では最大化問題だけを扱うため、以降の「上界」「下界」はすべてこの意味で使います。

この例では、単純な貪欲法が14点、後ほど説明するビームサーチが品物3と4を選んで20点、価格から計算した上界が25.492点になります。したがって、本当の最適値は20点以上25.492点以下です。

普通のDPではだめなのか

重さだけを制約にした1次元ナップサックなら、品物数を$N$、容量を$C$として、よく知られた$O(NC)$の動的計画法(DP)で解けます。$O(NC)$は、計算量がおおむね品物数$N$と容量$C$の積に比例して増えるという意味です。

重さと体積の2制約にそのまま拡張すると、計算量は

O(NC_{\mathrm{weight}}C_{\mathrm{volume}})

です。ここでは、品物を選ぶと消費され、使える量に上限があるものを資源と呼びます。この例では重さと体積が資源です。資源の種類数を$M$とすると、すべての容量の組合せを表に持つ素直なDPの状態数は、容量の積になります。

O\left(N\prod_{r=0}^{M-1}(C_r+1)\right)

容量が小さければDPを使うべきです。一方、容量が大きい、制約が5本や10本ある、重さは小さいが体積は大きい状態と、その逆の状態のどちらも捨てられず大量に残る、といった場合には、状態数を保持できません。この「一方の使用量では優れるが、別の使用量では劣るため単純には捨てられない状態」をPareto最適な状態と呼びます。

本記事の方法は、必ず最適解を返す厳密解法ではありません。ただし、計算量は容量の数値そのものではなく、各段階で残す途中状態の数と制約数に依存します。ビーム幅を$W$、上界の個数を$K$、1状態から試す行動数を$A$とすると、ビーム部分は大まかには、

O(NWAMK)

です。このナップサックでは$A=2$です。容量が1万から100万になっても、状態として持つ使用量は$M$個の整数のままです。

4品物の例で、入力から出力までを時系列に追う

ここから、この4品物・2資源の例を最初から最後まで追いながら、必要な概念も順番に説明します。価格調整とビームサーチを交互に行うのではなく、先に価格を決め、固定した価格でビームサーチを1回だけ行います。

1. 初期解を作る

先ほど説明した単純な貪欲法は品物1を選び、14点の完成解を作ります。この初期解は、

  • 価格更新の歩幅を決めるための「実際に14点は作れた」という基準
  • ビーム中に14点を超えられない状態を捨てるための基準
  • 最後に出力する候補

として使います。

2. 価格を調整し、3本を選んで固定する

まず、重さ1に$\lambda_w$、体積1に$\lambda_v$という0以上の仮の価格を付けます。品物$j$の、価格差し引き後の価値は、

p_j-\lambda_w a_{w,j}-\lambda_v a_{v,j}

です。容量制約を一時的に外すと、各品物を独立に考えられます。価格差し引き後の価値が正なら選び、負なら選ばない方法が、この制約を外した問題の最適解です。これは容量を守る先ほどの貪欲法とは別のものです。このように制約を外して簡単にした問題を緩和問題、制約の代わりに価格を付ける方法をラグランジュ緩和、仮の価格をラグランジュ乗数と呼びます。

価格$\boldsymbol{\lambda}=(\lambda_w,\lambda_v)$から、元の問題の上界を

U(\boldsymbol{\lambda})
=\lambda_w C_w+\lambda_v C_v
+\sum_{j=1}^{N}
\max(0,p_j-\lambda_w a_{w,j}-\lambda_v a_{v,j})

と計算します。最初の2項は、使える重量容量と体積容量全体の価格です。

この値が上界になる理由を確認します。容量を守る任意の選び方を$S$、その価値を$Z$とすると、

Z=\sum_{j\in S}p_j

です。各品物の価値を「使った資源の価格」と「価格差し引き後の価値」に分け、実際の使用量を容量まで増やしてよいことを使うと、

\begin{aligned}
Z
\le{}&\lambda_w C_w+\lambda_v C_v\\
&+\sum_{j\in S}
(p_j-\lambda_w a_{w,j}-\lambda_v a_{v,j})\\
\le{}&U(\boldsymbol{\lambda})
\end{aligned}

となります。最後の不等号では、残りの全品物についても価格差し引き後の価値が正なら足しています。したがって、どの容量内の選び方も$U(\boldsymbol{\lambda})$を超えず、元の最適値もこの値以下です。

価格更新の1反復は、次の順序です。

  1. 現在価格で、容量制約を外した問題を解く
  2. そこで使われた重量と体積を数える
  3. 容量との差から次の価格を決める

これを最大60回繰り返します。さまざまな価格を試して上界$U(\boldsymbol{\lambda})$を小さくする側の問題を双対問題と呼びます。制約を守る貪欲解は毎回作らず、最初に得た14点を、実際に作れるスコアの基準として使い続けます。

価格は$(0,0)$から始めます。60回すべてを並べると長いため、まず1回目と2回目を途中計算まで追い、その後の数回を表で確認します。

価格更新では、容量との差

\boldsymbol{g}
=(\text{使用重量}-10,\ \text{使用体積}-10)

をまず計算します。この$\boldsymbol{g}$が、価格を動かしたい生の方向です。ただし、そのまま使うと価格が激しく往復しやすいため、前回の更新方向も混ぜます。

過去方向を残す割合を$\alpha=0.35$とします。新しい生の方向へ掛ける0.65は別の設定値ではなく、$1-\alpha$です。数式では、

\boldsymbol{d}_{\mathrm{new}}
=0.35\,\boldsymbol{d}_{\mathrm{old}}+0.65\,\boldsymbol{g}

として得た$\boldsymbol{d}_{\mathrm{new}}$を実際の更新方向にします。
最初は$\boldsymbol{d}_{\mathrm{old}}=(0,0)$です。以下では、更新後の$\boldsymbol{d}_{\mathrm{new}}$を単に$\boldsymbol{d}$と書きます。

ここで、$\boldsymbol{\lambda}$と$\boldsymbol{d}$は別の量です。

  • $\boldsymbol{\lambda}$:現在の価格
  • $\boldsymbol{d}$:価格を動かす方向と、その相対的な強さ
  • $\Delta\boldsymbol{\lambda}$:今回の価格変化量

$\boldsymbol{d}$は名前どおり方向ですが、それ自体が価格変化量$\Delta\boldsymbol{\lambda}$なのではありません。歩幅$s$を掛けた

\Delta\boldsymbol{\lambda}=s\boldsymbol{d}

が価格変化量です。

更新後の価格が負になった場合は0へ戻します。容量が「上限以下」という制約では、上界の保証に0以上の価格が必要だからです。このように、値を許される範囲へ戻す処理を射影と呼びます。

また、$\boldsymbol{d}$を0から1の範囲へ収めたり、長さ1へ正規化したりする必要はありません。ここでいう「方向」は単位ベクトルという意味ではなく、各資源について価格を上げるか下げるかに加えて、過不足の相対的な大きさも持つベクトルです。そのため、成分は1より大きくなることも、価格を下げたいときには負になることもあります。

歩幅$s$は、$\boldsymbol{g}$と$\boldsymbol{d}$の内積に応じて調整します。仮に$\boldsymbol{d}$全体が同じ倍率で大きくなると、$s$はおおむねその倍率だけ小さくなるため、価格変化量$s\boldsymbol{d}$が必要以上に大きくなりません。

$\alpha$を大きくすると価格は滑らかに動きますが、現在の過不足への反応は遅くなります。小さくすると素早く反応しますが、価格が往復しやすくなります。0.35は本記事で固定した調整値です。

価格の更新は、

\boldsymbol{\lambda}_{\mathrm{next}}
=\boldsymbol{\lambda}+s\boldsymbol{d}
=\boldsymbol{\lambda}+\Delta\boldsymbol{\lambda}

です。歩幅$s$は、歩幅係数を$\eta$として、概ね次の式で決めます。

s
=\eta
\frac{U(\boldsymbol{\lambda})-\text{現在の最良解の価値}}
     {\boldsymbol{g}\cdot\boldsymbol{d}}

$\boldsymbol{g}\cdot\boldsymbol{d}$は、各成分の積を足したものです。上界と実際に作れた解の差が大きければ大きく動かし、方向の大きさに応じて割ることで動きすぎを抑えます。

歩幅係数$\eta$は1.15から始め、次の規則で更新します。

状況 歩幅係数$\eta$の更新 意味
これまでで最も小さい上界を得た $\eta\leftarrow\min(1.35,1.04\eta)$ 価格探索がうまく進んでいると判断し、今回の価格更新で使う係数を4%大きくする
現在価格の評価と次の価格への更新を3回行っても、上界の最小記録を更新できない $\eta\leftarrow\max(0.45,0.72\eta)$ 動きが大きすぎる可能性があるので、次の歩幅を28%小さくする
それ以外 変更しない 現在の大きさを維持する

ここでいう「改善」は、その反復の上界が直前の反復より小さいという意味ではなく、それまでの全反復で最も小さい上界を更新したという意味です。上界は小さいほど元の問題の最適値に近く、よい価格と判断できます。

停滞したときは歩幅係数を小さくするだけでなく、現在価格をそれまで最もよかった価格へ65%寄せ、更新方向を35%に弱めます。

この方向の弱化は、通常の平滑化を置き換えるものではなく、その反復で通常更新した後に追加します。再中心化では現在価格を35%残して最良価格を65%混ぜ、更新方向も35%残します。この35%と、通常更新で過去方向を残す割合$\alpha=0.35$は、同じ数値ですが独立した調整値です。

停滞処理を行った回を新しい起点とし、そこから「現在価格を評価して次の価格へ更新する」という一巡をもう一度数え始めます。その一巡を3回行っても最小記録を更新できなければ、再び停滞処理を行います。途中で最小記録を更新できた場合は、その改善した回から数え直します。たとえば7回目に停滞処理を行い、8、9、10回目のどこでも最小記録を更新できなければ、10回目にも停滞処理を行います。

歩幅係数$\eta$と、最終的な歩幅$s$は別の量です。$s$は上の式で毎回計算した後、設定した最小値以上、最大値以下に制限します。この例では最小値を0.01とし、最大値を価格調整の序盤の10.0から終盤の1.0まで直線的に小さくします。

1回目では、価格は$\boldsymbol{\lambda}=(0,0)$です。価格が0なので、4品物の価格差し引き後の価値はすべて正です。容量制約を外した緩和問題では、4品物をすべて選びます。

U(\boldsymbol{\lambda})=14+14+10+10=48
\qquad(\boldsymbol{\lambda}=(0,0))
項目 計算方法
選ぶ品物 ${1,2,3,4}$ 価格が0なので、価格差し引き後の価値$14,14,10,10$がすべて正
使用量 $(21,21)$ $(10+1+5+5,\ 1+10+5+5)$
容量 $(10,10)$ 入力で与えられた上限
生の更新方向$\boldsymbol{g}$ $(11,11)$ $(21,21)-(10,10)$
過去方向を混ぜた$\boldsymbol{d}$ $(7.15,7.15)$ $0.35(0,0)+0.65(11,11)$
上界$U(\boldsymbol{\lambda})$ 48.000 $0+14+14+10+10$
実際に作れた初期解 14 単純な貪欲法で品物1を選んだ価値
上界と初期解の差 34.000 $48-14$
歩幅係数$\eta$ 1.196 $\min(1.35,\ 1.15\times1.04)$
歩幅$s$ 0.258512 $\displaystyle 1.196\frac{48-14}{11\times7.15+11\times7.15}$
更新後の価格 $(1.848364,1.848364)$ $(0,0)+0.258512397(7.15,7.15)\simeq(1.848364,1.848364)$

容量を11ずつ超えたため、重さと体積の価格がどちらも大きく上がりました。

2回目では、価格は$\boldsymbol{\lambda}=(1.848364,1.848364)$です。今度は価格が高くなりすぎました。たとえば品物1の価格差し引き後の価値は、

14-1.848364\times10-1.848364\times1<0

です。他の品物もすべて負になるため、緩和問題は何も選びません。品物から来る$\max(0,\cdot)$の項はすべて0になり、上界には容量全体の価格だけが残ります。

U(\boldsymbol{\lambda})
=10\times1.848364+10\times1.848364
\simeq36.967
\qquad(\boldsymbol{\lambda}=(1.848364,1.848364))
項目 計算方法
選ぶ品物 なし $14-10\times1.848364-1\times1.848364<0$かつ$10-5\times1.848364-5\times1.848364<0$
使用量 $(0,0)$ 何も選ばないため
容量 $(10,10)$ 入力で与えられた上限
生の更新方向$\boldsymbol{g}$ $(-10,-10)$ $(0,0)-(10,10)$
過去方向を混ぜた$\boldsymbol{d}$ $(-3.9975,-3.9975)$ $0.35(7.15,7.15)+0.65(-10,-10)$
上界$U(\boldsymbol{\lambda})$ 36.967 $10\times1.848364+10\times1.848364$
上界と初期解の差 22.967 $36.967273-14$
歩幅係数$\eta$ 1.24384 $\min(1.35,\ 1.196\times1.04)$
歩幅$s$ 0.357318 $\displaystyle 1.24384\frac{36.967273-14}{(-10)(-3.9975)+(-10)(-3.9975)}$
更新後の価格 $(0.419983,0.419983)$ $(1.848364,1.848364)+0.35731848(-3.9975,-3.9975)\simeq(0.419983,0.419983)$

今回は資源をまったく使わず10ずつ余らせたため、価格を下げます。

最初の6回をまとめると、次のようになります。「更新前価格」で緩和問題を解き、右端の「更新後価格」を次の反復で使います。

反復 更新前価格$\boldsymbol{\lambda}$ 緩和問題で選ぶ品物 使用量 生の方向$\boldsymbol{g}$ 上界$U(\boldsymbol{\lambda})$ 更新後価格
1 $(0.000,0.000)$ ${1,2,3,4}$ $(21,21)$ $(11,11)$ 48.000 $(1.848,1.848)$
2 $(1.848,1.848)$ なし $(0,0)$ $(-10,-10)$ 36.967 $(0.420,0.420)$
3 $(0.420,0.420)$ ${1,2,3,4}$ $(21,21)$ $(11,11)$ 38.760 $(1.820,1.820)$
4 $(1.820,1.820)$ なし $(0,0)$ $(-10,-10)$ 36.398 $(0.371,0.371)$
5 $(0.371,0.371)$ ${1,2,3,4}$ $(21,21)$ $(11,11)$ 39.833 $(1.890,1.890)$
6 $(1.890,1.890)$ なし $(0,0)$ $(-10,-10)$ 37.804 $(0.351,0.351)$

この6反復で使われる歩幅係数$\eta$は、順に1.196、1.24384、1.24384、1.293594、1.293594、1.293594です。1回目は最初の上界なので最小記録になり、2回目は48.000から36.967へ、4回目は36.967から36.398へ最小記録を更新したため、それぞれ1.04倍されます。3、5、6回目は最小記録を更新していないため変わりません。4回目の改善後、6回目の終了時点ではまだ停滞が2反復なので、0.72倍する条件には達していません。

このように、序盤は「安すぎて全部選ぶ」と「高すぎて何も選ばない」の間を往復します。そこでライブラリは、過去の更新方向を混ぜ、改善が止まったらそれまで最も上界が小さかった価格の近くへ戻し、歩幅も調整します。各反復の価格と上界は、3本を後で選ぶためのスナップショットとして保存します。

60回の更新後、保存したスナップショットから3本を選びます。まず、初期状態の上界が最も小さい価格を1本目にします。残りは上界が小さい上位14本を候補とし、「すでに選んだ価格から遠く、上界もあまり悪化しない」ものを追加します。

候補$q$から選択済み価格への最短距離を$d_q$、最良上界からの悪化量を$\ell_q$とし、

\frac{d_q}{1+\ell_q/5000}

が最も大きい候補を順に選びます。距離は資源ごとの価格の大きさで正規化してから測るため、単位の違う資源も比較できます。

初期状態で最良の価格が、すべての途中状態でも最良とは限りません。使った資源量と残り品物が変わると、途中状態の上界を最も小さくする価格も変わるためです。そこで、ビームサーチでは3本それぞれから上界を計算し、その最小値を使います。どの価格から得た値も上界なので、その最小値も安全な上界です。

この小入力で実際に選ばれた3本は次の通りです。数値は小数第6位まで示しています。

番号 重さの価格$\lambda_w$ 体積の価格$\lambda_v$ 初期状態の上界$U(\boldsymbol{\lambda})$ 選ばれた理由
0 1.253752 1.253752 25.492496 上界が最小だった1本目
1 0.928542 0.928542 約27.572 1本目から離れた良質な候補
2 1.078089 1.078089 約25.844 既選択の価格と上界の両方を考えて追加

この問題は重さと体積について対称なので、3本とも$\lambda_w=\lambda_v$になっています。3本を選んだ時点で価格調整は終わり、以後は固定します。この後のビームサーチ中に価格は変えません。

3. ビームサーチを始める前に上界を前計算する

資源を$r=0,1,\ldots,M-1$と番号付けし、資源$r$の容量を$C_r$、品物$j$が使う量を$a_{r,j}$、価格を$\lambda_r$とします。ビームサーチの途中では、すでに選んだ品物の価値を$V$、資源$r$の使用量を$u_r$、まだ処理していない品物の集合を$R$とします。価格$\boldsymbol{\lambda}$に対する途中状態の上界は、問題全体の$U(\boldsymbol{\lambda})$を確定済み部分と残り部分に分けた、

B(\boldsymbol{\lambda})
=V
+\sum_{r=0}^{M-1}\lambda_r(C_r-u_r)
+\sum_{j\in R}
\max\left(0,p_j-\sum_{r=0}^{M-1}\lambda_r a_{r,j}\right)

です。「現在の価値」「残り容量の価格」「残り品物を制約なしで選んだ価値」の和です。初期状態では$V=0$、$u_r=0$、$R$が全品物なので、問題全体の上界$U(\boldsymbol{\lambda})$と一致します。途中状態から容量を守って完成させたどの解も、この値を超えません。

各価格について、まだ見ていない品物から得られる「価格差し引き後の正の価値」の合計を、後ろから累積して表にします。これに現在価値と残り容量の価格を足せば、各候補の上界をすぐ計算できます。

まず、価格番号$k$と品物$j$の組について、価格差し引き後の正の価値

q_{k,j}
=\max\left(0,
p_j
-\underbrace{\lambda_{k,w}a_{w,j}}_{\text{重量価格$\times$重量}}
-\underbrace{\lambda_{k,v}a_{v,j}}_{\text{体積価格$\times$体積}}
\right)

を計算します。$k$は、先ほど選んだ3本の価格の番号です。この例では重さと体積の価格がたまたま同じ数値ですが、どちらの資源に対する項か分かるよう、以下では「重量価格$\times$重量」と「体積価格$\times$体積」を分けて書きます。

価格番号$k$ 価格$(\lambda_{k,w},\lambda_{k,v})$ 品物1の$q_{k,1}$ 品物2の$q_{k,2}$ 品物3の$q_{k,3}$ 品物4の$q_{k,4}$
0 $(1.253752,1.253752)$ $\max(0,14-1.253752\times10-1.253752\times1)=0.208728$ $\max(0,14-1.253752\times1-1.253752\times10)=0.208728$ $\max(0,10-1.253752\times5-1.253752\times5)=0$ $\max(0,10-1.253752\times5-1.253752\times5)=0$
1 $(0.928542,0.928542)$ $14-0.928542\times10-0.928542\times1=3.786038$ $14-0.928542\times1-0.928542\times10=3.786038$ $10-0.928542\times5-0.928542\times5=0.714580$ $10-0.928542\times5-0.928542\times5=0.714580$
2 $(1.078089,1.078089)$ $\max(0,14-1.078089\times10-1.078089\times1)=2.141021$ $\max(0,14-1.078089\times1-1.078089\times10)=2.141021$ $\max(0,10-1.078089\times5-1.078089\times5)=0$ $\max(0,10-1.078089\times5-1.078089\times5)=0$

ビームサーチでは、処理済みの品物をこの合計から外していきます。毎回足し直さずに済むよう、まだ処理していない品物の$q_{k,j}$の合計を前計算します。

ビームサーチの位置 まだ処理していない品物 価格0での残り合計 価格1での残り合計 価格2での残り合計
開始前 1、2、3、4 0.417456 9.001236 4.282042
品物1の処理後 2、3、4 0.208728 5.215198 2.141021
品物2の処理後 3、4 0 1.429160 0
品物3の処理後 4 0 0.714580 0
品物4の処理後 なし 0 0 0

ある途中状態について、価格$k$から計算する上界を$B_k$と書くと、使う式は次の一つです。

B_k
=\text{現在価値}
+\underbrace{\lambda_{k,w}(10-\text{使用重量})}_{\text{残り重量の価格}}
+\underbrace{\lambda_{k,v}(10-\text{使用体積})}_{\text{残り体積の価格}}
+\text{価格$k$での残り合計}

最後に3本から得た$B_0,B_1,B_2$のうち、最も小さいものを、その状態の上界として表に載せています。

たとえば初期状態は現在価値0、使用量$(0,0)$で、全品物が残っています。3本それぞれでは、

\begin{aligned}
B_0&=0+1.253752\times10+1.253752\times10+0.417456=25.492496,\\
B_1&=0+0.928542\times10+0.928542\times10+9.001236=27.572076,\\
B_2&=0+1.078089\times10+1.078089\times10+4.282042=25.843822
\end{aligned}

となります。最も小さい$B_0=25.492496$を選ぶため、初期状態の表では上界を25.492と表示しています。前計算はこの式を高速に評価するためのもので、別の上界を計算しているわけではありません。

この例では、初期解の14点を直後のビームサーチにおける枝刈りの基準にします。

途中状態の上界が14点を下回るなら、そこからどのように完成させても初期解を超えられないため、その状態は安全に捨てられます。この処理を数理最適化ではincumbent cutoffと呼びます。一方、上界が高くてもビーム幅$W$から漏れて捨てられる状態はあり、こちらは計算量を抑えるためのヒューリスティックです。

以下の表では、状態を次の形で表します。

{選んだ品物}; 現在価値; 使用量(重さ, 体積); 3価格から得た上界の最小値

初期状態は、

選んだ品物 現在価値 使用量 上界
$\varnothing$ 0 $(0,0)$ 25.492

です。この例ではビーム幅を3にします。候補が4状態以上になれば、容量を守り上界が初期解以上でも、上界が大きい上位3状態だけを残します。

一つの候補については、3価格から得た上界の最小値を使い、できるだけ厳しく見積もります。一方、異なる候補どうしを並べるときは、その上界が大きい候補を「将来高得点になる可能性がある」と見て優先します。実際の大規模問題では、この順に上位$W$状態だけを残します。

計算量を抑えるため、まず1本目の価格だけで全候補を評価して上位$2W$個へ絞り、その候補だけを残り2本でも評価します。最後に3本の最小値で並べ直し、上位$W$個を残します。

4. 品物1を処理する

初期状態から、品物1を「選ばない」「選ぶ」の2候補を作ります。

「選ばない」場合は現在価値と使用量がそのままで、残り品物が2、3、4になります。

\begin{aligned}
B_0&=0+1.253752\times10+1.253752\times10+0.208728=25.283768,\\
B_1&=0+0.928542\times10+0.928542\times10+5.215198=23.786038,\\
B_2&=0+1.078089\times10+1.078089\times10+2.141021=23.702801
\end{aligned}

したがって、この候補の上界は3本の最小値23.702801、表の表示では23.703です。

「選ぶ」場合は品物1の価値14が現在価値へ入り、使用量が$(10,1)$になります。品物1は処理済みなので、こちらでも残り品物は2、3、4です。

\begin{aligned}
B_0&=14+1.253752\times0+1.253752\times9+0.208728=25.492496,\\
B_1&=14+0.928542\times0+0.928542\times9+5.215198=27.572076,\\
B_2&=14+1.078089\times0+1.078089\times9+2.141021=25.843822
\end{aligned}

したがって、この候補の上界は3本の最小値25.492496、表の表示では25.492です。品物1を選んだ価値14は最初の項へ移り、使った重量10と体積1は残り容量をそれぞれ0と9へ減らしています。

品物2以降の表もまったく同じ計算です。「選ばない」場合は現在価値と使用量を変えず、残り合計からその品物を外します。「選ぶ」場合は品物の価値を現在価値へ足し、重さと体積を使用量へ足したうえで、その品物を残り合計から外します。その状態を$B_0,B_1,B_2$へ代入し、最小値を表の上界にします。

親状態 行動 選んだ品物 現在価値 使用量 上界 判定
$\varnothing$ 選ばない $\varnothing$ 0 $(0,0)$ 23.703 残す
$\varnothing$ 選ぶ ${1}$ 14 $(10,1)$ 25.492 残す

上界が大きい順では${1}$、$\varnothing$です。どちらも初期解の14点以上へ到達できる可能性があるため残ります。

5. 品物2を処理する

残った2状態から、それぞれ2候補ずつ、合計4候補を作ります。

親状態 行動 選んだ品物 現在価値 使用量 上界 判定
$\varnothing$ 選ばない $\varnothing$ 0 $(0,0)$ 20.000 残す
$\varnothing$ 選ぶ ${2}$ 14 $(1,10)$ 23.703 残す
${1}$ 選ばない ${1}$ 14 $(10,1)$ 23.703 残す
${1}$ 選ぶ ${1,2}$ 28 $(11,11)$ 容量違反なので作らない

品物1と2を両方選ぶ候補は、重さも体積も上限10を超えるため除外します。次の段階へ進むのは$\varnothing,{1},{2}$の3状態です。

6. 品物3を処理する

親状態 行動 選んだ品物 現在価値 使用量 上界 判定
$\varnothing$ 選ばない $\varnothing$ 0 $(0,0)$ 19.285 ビーム幅で捨てる
$\varnothing$ 選ぶ ${3}$ 10 $(5,5)$ 20.000 残す
${1}$ 選ばない ${1}$ 14 $(10,1)$ 23.071 残す
${1}$ 選ぶ ${1,3}$ 24 $(15,6)$ 重量違反なので作らない
${2}$ 選ばない ${2}$ 14 $(1,10)$ 23.071 残す
${2}$ 選ぶ ${2,3}$ 24 $(6,15)$ 体積違反なので作らない

容量を守る4候補ができましたが、ビーム幅は3です。上界が大きい${1},{2},{3}$を残し、上界19.285の$\varnothing$を捨てます。${3}$は現在10点しかありませんが、品物4も選べば20点になるため、上界20として残ります。

7. 品物4を処理する

親状態 行動 選んだ品物 現在価値 使用量 上界 判定
${1}$ 選ばない ${1}$ 14 $(10,1)$ 22.357 残す
${1}$ 選ぶ ${1,4}$ 24 $(15,6)$ 重量違反なので作らない
${2}$ 選ばない ${2}$ 14 $(1,10)$ 22.357 残す
${2}$ 選ぶ ${2,4}$ 24 $(6,15)$ 体積違反なので作らない
${3}$ 選ばない ${3}$ 10 $(5,5)$ 19.285 ビーム幅で捨てる
${3}$ 選ぶ ${3,4}$ 20 $(10,10)$ 20.000 残す

容量違反を除くと4候補ができます。上界が大きい${1},{2},{3,4}$を残し、上界19.285の${3}$をビーム幅によって捨てます。

8. 完成状態を本来の価値で比較する

すべての品物について判断した後も、上界と本来の価値は一致しないことがあります。たとえば${1}$の本来の価値は14ですが、未使用の体積容量に価格が残るため、表の上界は22.357です。

そこで、最後は上界ではなく本来の価値で比較します。

完成状態 本来の価値 最終候補としての結果
${1}$ 14 選ばない
${2}$ 14 選ばない
${3,4}$ 20 ビーム解として選ぶ

最後に初期解14点とビーム解20点を本来の価値で比較し、ビーム解${3,4}$を最終出力にします。価格調整で得た緩和解や上界そのものを回答として出力することはありません。

ライブラリと問題固有部分の境界

ここからは、ここまで説明した処理をライブラリへどう分けたかを説明します。問題ごとに用意するのは、「価格を固定した小問題の解き方」と「解をどの順に構築するか」の二つです。ライブラリはそれらを受け取り、価格調整とビームサーチを行います。

ライブラリが担当するもの

問題に依存しない共通部分は、次を担当します。

  1. 決められた回数だけ価格を更新する
  2. 上界と現在の最良解の差から価格の歩幅を決める(Polyak型ステップ)
  3. 過去の更新方向も混ぜ、価格の往復を抑える(価格安定化)
  4. 改善が止まったら、これまで最もよかった価格の近くへ戻す(再中心化)
  5. 上界と、それを計算した価格ベクトルを保存する
  6. 上界が小さく、互いに似すぎていない価格を選ぶ
  7. 複数の上界を使ってビームの途中状態を並べる
  8. 現在の最良解を超えられない状態の削除
  9. 同じ状態ハッシュを持つ候補を一つにまとめる(ハッシュによる重複除去)
  10. 「どの行動を選んだか」という履歴から最終回答を復元する

時間打ち切りはありません。価格を更新する回数とビーム幅を固定すれば、探索量は入力と設定値だけで決まります。

価格調整について問題側が決めるもの

価格を固定したときに、元の問題から分かれてできる小問題を定義します。たとえばナップサックでは各品物を選ぶかどうか、Introでは各コンテスト種類をどの日に開催するかが小問題です。

  • 価格の初期値
  • 制約を守る最初の解
  • 容量と価格だけで決まる定数項
  • 価格を固定した各小問題の最良値
  • 制約の過不足から価格を動かす方向を求める方法
  • 価格が許容範囲を出たときに戻す方法

ナップサックでは1品物が1つの独立部分問題です。Introでは1コンテスト種類の全365日計画が1つの独立部分問題です。

ビームサーチについて問題側が決めるもの

解をどの順番で作るかを定義します。

  • 各段階で試す行動
  • その行動が制約を守るか
  • 行動による確定済みスコアと状態の変化
  • 同じ状態へ到達した候補を見分ける識別値
  • 価格差し引き後の未来値と、その差分更新
  • 未確定部分について上界へ足し戻す価格
  • 完成状態の正確な採点と回答形式への変換

サンプルコードでは、これらをProblemSearchModelに分けて実装しています。必要なインターフェースはC++20のconceptで検査するため、関数の不足や型の不一致はコンパイル時に検出されます。

※本記事の方法はラグランジュ緩和、双対価格の調整、現在の最良解による枝刈りなど、数理最適化で使われる考え方を取り入れています。この上界を利用する別の探索として、上界の高い状態から展開するbest-first分枝限定も試しましたが、探索が広がりすぎて初期解を改善できなかったため、採用しませんでした。

同じライブラリでIntroを解く

ここからは、同じ共通部分をIntroへ適用します。問題固有部分だけを交換し、価格調整とビームサーチの本体は変更しません。

Introの難しさを分けて考える

Introでは、365日それぞれに26種類のうち1種類のコンテストを割り当てます。

  • その日に開催した種類の満足度を得る
  • 各種類を長く開催しないほど、毎日の不満ペナルティが増える
  • 1日に開催できる種類はちょうど1つ

一見すると26種類が複雑に影響し合っています。しかし、図のように、ある1種類の開催日だけに注目すると、その種類の寄与は独立に計算できます。

intro_single_contest_independence.png

コンテスト種類の番号を$t$、種類$t$の不満係数を$c_t$とします。種類$t$を開催しない期間が$n$日続いたときのペナルティは、

c_t(1+2+\cdots+n)
=c_t\frac{n(n+1)}2

です。以下、この三角数を$\operatorname{tri}(n)=n(n+1)/2$と書きます。

種類ごとに独立でない原因は、ただ一つです。

どの日も、26種類のうちちょうど1種類だけを開催する。

そこで、この制約を一時的に外します。

各日に場所代を付ける

日$d$に価格$\lambda_d$を付けます。26種類はそれぞれ独立に、自分をいつ開催すると最も得かを考えます。

種類$t$の開催日集合を$S_t$、種類$t$だけの寄与を$F_t(S_t)$と書くと、小問題は、

\max_{S_t}
\left(
F_t(S_t)-\sum_{d\in S_t}\lambda_d
\right)

です。問題全体の上界は、全日の場所代を足し戻した、

U(\boldsymbol{\lambda})
=\sum_{d=0}^{D-1}\lambda_d
+\sum_{t=0}^{25}
\max_{S_t}
\left(
F_t(S_t)-\sum_{d\in S_t}\lambda_d
\right)

です。正しい日程では各日をちょうど1種類が使うため、引いた場所代と第1項が相殺します。一方、制約を外した小問題では各種類が都合のよい日程を選べるため、この値は正しい日程のスコア以上、つまり上界になります。「1日ちょうど1種類」は等式制約なので、$\lambda_d$は正でも負でも構いません。

1種類だけならDPで解ける

価格を固定した1種類の問題はDPで解けます。

日を$0,1,\ldots,D-1$の0始まりで表します。$dp[d]$を「最後にその種類を開催した日が$d$であるときの、日$d$までの価格差し引き後の最大値」とします。$dp[d]$には日$d$より後の不満はまだ含めません。式中の記号は次の意味です。

  • $d$:今回その種類を開催する日
  • $l$:その一つ前に開催した日
  • $s_{d,t}$:日$d$に種類$t$を開催したときの満足度
  • $c_t$:種類$t$を開催しないことで増える、1日あたりの不満係数
  • $\operatorname{tri}(n)=n(n+1)/2$:$1+2+\cdots+n$を返す三角数
  • $\lambda_d$:日$d$を1種類が使うことへの仮の価格

まず、4日間で1種類だけを考える小さな例を順に追います。不満係数を$c_t=1$とし、4日間の(満足度、価格)を順に$(8,1),(1,4),(9,2),(1,4)$とします。

各日について、「その日が最初の開催日である場合」と、「それ以前の各日を直前の開催日とする場合」を比べます。開催しなかった日が間に$n$日あれば、不満$\operatorname{tri}(n)$を引きます。

今回の開催日$d$ その日が最初の開催 直前の開催日$l$ごとの候補 $dp[d]$
0 $8-1-\operatorname{tri}(0)=7$ なし 7
1 $1-4-\operatorname{tri}(1)=-4$ $l=0:\ dp[0]+1-4-\operatorname{tri}(0)=4$ 4
2 $9-2-\operatorname{tri}(2)=4$ $l=0:\ dp[0]+9-2-\operatorname{tri}(1)=13$
$l=1:\ dp[1]+9-2-\operatorname{tri}(0)=11$
13
3 $1-4-\operatorname{tri}(3)=-9$ $l=0:\ dp[0]+1-4-\operatorname{tri}(2)=1$
$l=1:\ dp[1]+1-4-\operatorname{tri}(1)=0$
$l=2:\ dp[2]+1-4-\operatorname{tri}(0)=10$
10

たとえば$dp[2]=13$は、日0と日2に開催する場合です。開催による満足度から価格を引いた値は$(8-1)+(9-2)$、間の日1に開催しない不満が1なので、$(8-1)+(9-2)-1=13$になります。

一般に、直前の開催日を$l$とすると、同じ計算は次の式で表せます。

\begin{aligned}
dp[d]
={}&\max\left(
s_{d,t}-\lambda_d-c_t\operatorname{tri}(d),\right.\\
&\left.\qquad
\max_{l<d}\left\{
dp[l]+s_{d,t}-\lambda_d
-c_t\operatorname{tri}(d-l-1)
\right\}
\right)\\
={}&s_{d,t}-\lambda_d
+\max\left(
-c_t\operatorname{tri}(d),
\max_{l<d}\left\{
dp[l]-c_t\operatorname{tri}(d-l-1)
\right\}
\right).
\end{aligned}

1本目の等式が表の計算に対応し、2本目では両方の候補に共通する$s_{d,t}-\lambda_d$をくくり出しています。外側の$-c_t\operatorname{tri}(d)$は、日$d$が最初の開催日で、それより前の$d$日間に一度も開催しなかった場合です。内側は、直前の開催日$l$から今回の開催日$d$までの間に、$d-l-1$日間開催しなかった場合です。

この時点では、最後の開催日より後の不満がまだ入っていません。4日間の例では、最後の開催日ごとの完成値は次のようになります。

最後の開催日 それ以降の不満を引く計算 完成値
一度も開催しない $-\operatorname{tri}(4)$ -10
0 $dp[0]-\operatorname{tri}(3)=7-6$ 1
1 $dp[1]-\operatorname{tri}(2)=4-3$ 1
2 $dp[2]-\operatorname{tri}(1)=13-1$ 12
3 $dp[3]-\operatorname{tri}(0)=10-0$ 10

最大は12なので、この種類は日0と日2に開催するのが最善です。これで、最初の開催前、開催日の間、最後の開催後の不満がすべて入りました。

一般には、日$D-1$まで一度も開催しない場合と、最後に日$d$で開催する場合を比較します。種類$t$の小問題の最良値$V_t$は、

V_t
=\max\left(
-c_t\operatorname{tri}(D),
\max_{0\le d<D}
\left\{
dp[d]-c_t\operatorname{tri}(D-d-1)
\right\}
\right)

です。

全日数を$D$とすると、1種類あたり$O(D^2)$、26種類でも$O(26D^2)$です。Introでは$D=365$なので、価格を変えながら何度も解けます。

人気日の価格を上げる

26種類を独立に解いた後、各日が何種類から選ばれたかを数えます。

日$d$を開催日に選んだ種類数を$q_d$とすると、価格を動かす生の方向は$q_d-1$です。本来は各日をちょうど1種類が使うため、目標値の1を引いています。

たとえば3種類が同じ日を選んだなら$3-1=2$なので、その日の価格を上げます。どの種類も選ばなければ$0-1=-1$なので、価格を下げます。

  • 2種類以上から選ばれた日:価格を上げる
  • 1種類から選ばれた日:そのまま
  • どの種類からも選ばれなかった日:価格を下げる

ナップサックでは更新後の価格を0以上へ射影しました。一方、Introの「ちょうど1種類」は等式制約なので、価格は負でも構いません。そのため射影は必要ありません。

Intro固有なのは$q_d-1$という過不足の計算です。過去方向を混ぜる処理、歩幅、安定化、保存する価格の選択はナップサックと共通です。

Introのビーム状態

Introでは、日を前から1日ずつ決めます。

1段階を1日に対応させ、その日に開催する種類を一つ選びます。途中状態には、各種類を最後に開催した日と、ここまでに確定した満足度・不満を持たせます。

各日では、評価の高い6〜10種類を、その日に開催する種類の候補として選びます。候補ごとに、その種類を開催した子状態を一つずつ作ります。

親状態
├── 候補Aをその日に開催した子状態
├── 候補Bをその日に開催した子状態
└── 候補Cをその日に開催した子状態

「その日は何も開催しない」という選択肢はありません。そのため、探索の構造だけで「各日にちょうど1種類」を満たします。Introには開催可能日や開催回数に関する追加制約もないので、候補に入った種類はすべて選べます。ナップサックのような容量検査は不要です。

残り期間の上界を作る

途中状態では、各種類について「最後に開催した日」が分かっています。価格を固定し、次の日から365日目までの最良開催計画を種類別に解き、未確定日の場所代を足し戻せば、残り期間の上界になります。

この値を毎回DPで計算すると遅いため、未来の値を表にして先に計算しておきます。表は次の4項目を指定すると値が一つ決まる形です。

保存した価格の番号 × 種類 × 次の日 × 最終開催日

26種類の未来値を合計し、未確定日の価格を足し戻します。途中状態の上界は、

確定済みスコア + 種類別未来値の合計 + 未確定日の価格

です。複数の価格を使う場合は、価格ごとにこの値を計算し、その最小値を候補の上界にします。

ある日で種類$t$を開催したときは、26種類全部を再計算する必要はありません。親状態の未来値から種類$t$の古い未来値を引き、新しい最終開催日に対する未来値を足します。

更新後の未来値合計 = 更新前の未来値合計 - 古い未来値 + 新しい未来値

26種類分を毎回計算し直さず、今回開催した1種類の古い未来値を引き、新しい未来値を足して差分更新します。種類別未来値の合計だけでは上界にならず、確定済みスコアと未確定日の価格も必要です。

各日の候補を6~10種類へ絞る

26種類すべてをビーム展開すると重いため、最良の価格ベクトルを使い、日$d$に種類$t$を開催する有望さ$G_{d,t}$を計算します。

これは「この日に必ず$t$を開催する場合の最良値-何も固定しない場合の最良値」です。0以下の値になり、0に近い種類ほど有望です。

各日について$G_{d,t}$のよい順に並べ、最良値から5,000点以内にある種類数を数え、その個数を6以上10以下に収めます。この候補生成、種類別DP、三角数ペナルティはIntro固有です。

完成状態は本来のスコアで選ぶ

上界は、完成状態でも本来の目的値に一致するとは限りません。ナップサックでは、使わなかった容量の価格が上界に残ります。

そのため、最後は上界ではなく本来のスコアで比較します。途中までのスコアには、各種類を最後に開催してから最終日までの不満がまだ入っていません。種類$t$の最終開催日を$l_t$とすると、最後に$c_t\operatorname{tri}(D-l_t-1)$を引きます。一度も開催していない種類では$l_t=-1$とすれば、全期間の不満を同じ式で計算できます。

このように、上界の作り方と本来の採点を分離することで、不等式制約のナップサックと等式制約のIntroを同じ探索本体で扱えます。

Introでの設定

今回のIntroでは次の固定設定を使います。

項目 意味
価格更新 使用量の過不足を見て価格を調整する回数の上限 60回
保存して使う価格 上界の計算に同時に使う価格ベクトル数 3本
価格候補プール 60回の更新中から、保存候補として一時的に残す本数 14本
ビーム幅 各日を決めた後に残す途中状態数 14,500
1本目の上界による事前選抜 まず計算の軽い上界1本で候補をこの数まで絞る ビーム幅の2倍
重複除去表 同じ状態をまとめるためのハッシュ表の要素数 32,768
時間打ち切り 実行時間を見て途中終了する処理 なし
焼きなまし・局所探索 完成解を少しずつ変更する後処理 なし
DP後処理 完成解の一部をDPで作り直す後処理 なし

サンプルコードの使い方

ライブラリ、多次元0/1ナップサックの使用例、Introの使用例はGitHubリポジトリで公開しています。

リポジトリを取得し、多次元0/1ナップサックの例を実行するには、次のようにします。

git clone https://github.com/thun-c/intro-heuristics-lagrangian-bound-beam.git
cd intro-heuristics-lagrangian-bound-beam
g++ -std=gnu++23 -O2 -Isrc src/01_multidimensional_knapsack.cpp -o a.out
./a.out < examples/knapsack_small.txt

Intro版も同様です。

g++ -std=gnu++23 -O2 -Isrc src/02_intro_heuristics.cpp -o a.out
./a.out < input.txt > output.txt

AtCoderへ提出するために1ファイルへまとめる場合は、combiner.shを使います。

./combiner.sh src/02_intro_heuristics.cpp
g++ -std=gnu++23 -O2 combined/combined_02_intro_heuristics.cpp -o a.out

結合したコードはcombined/combined_02_intro_heuristics.cppへ生成されます。

このライブラリが向く問題、向かない問題

向いているのは次のような問題です。

  • ある制約を外すと、独立な小問題へ分かれる
  • 小問題を価格固定で高速に解ける
  • 制約の過不足を数えられる
  • 解を段階的に構築できる
  • 途中状態からの安全な上界、または下界を計算できる

一方、次の場合には適していません。

  • 緩和しても小問題が難しい
  • 上界が非常に緩く、ビーム順位に情報を与えない
  • 状態や行動を固定順序で構築しにくい
  • 小さなDPや通常の全探索で厳密に解ける
  • 問題側で、計算値が本当に上界または下界になることを数式で保証できない

特に重要なのは、「確定済みスコア+価格差し引き後の未来値+未確定部分で足し戻す価格」が、本当に最終スコアを超えない値になっていることです。保証のない予測値を枝の削除に使うと、本当は最良解へ到達できる状態まで誤って削除する可能性があります。プログラムは関数の型を検査できても、上界を導く数式の正しさまでは検査できません。

まとめ

本記事の方法では、問題構造から安全な未来上界を作り、その値でビームサーチを案内します。

  1. 結合制約に価格を付ける
  2. 制約を外して独立な小問題へ分ける
  3. 小問題から上界と価格更新方向を得る
  4. 性格の異なる複数の価格を残す
  5. 複数上界でビームを順位付けする
  6. 現在の最良解を超えられない途中状態を安全に捨てる

ナップサックでは、品物を独立にするために容量制約へ価格を付けました。Introでは、26種類を独立にするために各日の「ちょうど1種類」という制約へ価格を付けました。

問題固有の小問題と、状態・行動・未来上界を実装すれば、価格更新、安定化、複数上界、ビーム、枝刈り、重複除去には同じ共通ライブラリを利用できます。

筆者自身も、他のAHCでどれぐらいこの方針が役に立つかは確認できていません。しかし、ナップサックの例にあるように、うまく同じ形として扱えれば汎用性がありそうでした。読者の皆さんがこの方式をうまく活用し、AHCでいい順位をとってくれることを期待します。

1
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
1
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?