補遺: 量子アニーリング応用事例とQUBO設計

量子アニーリングを現場へ持っていくとき、最初に必要になるのは「量子らしさ」の議論ではありません。現場の言葉で書かれた制約、優先順位、運用ルール、評価指標を、0-1変数の二次式へ翻訳する作業です。ここがうまくいくと、量子アニーラ、シミュレーテッドアニーリング、量子モンテカルロ、デジタルアニーラ、GPU上のIsingソルバなど、いろいろなサンプラを同じ土俵で比べられます。

QUBOは、万能の魔法の箱ではありません。むしろ、現場の問題をいったん裸にするための形式です。何を選ぶのでしょうか。何を同時に選んではいけないのでしょうか。少し悪い解でも許せるのでしょうか。絶対に守る制約はどれか。こうした判断が、係数行列の中に入ります。量子アニーリングの応用研究の面白さは、物理の道具を使いながら、実はこの翻訳の深さで勝負しているところにあります。

1. QUBOという翻訳形式

0-1変数を $\mathit{x}=(x_1,\dots ,x_N)$、$x_i\in \{0,1\}$ とします。QUBOは、目的関数を

$$E(\mathit{x})={\mathit{x}}^{\mathsf{T}}Q\mathit{x}=\sum _i{Q}_{ii}{x}_i+\sum _{i<j}{Q}_{ij}{x}_i{x}_j$$

の形で書く最適化問題です。$Q_{ii}$ は「その選択肢を単独で選んだときの得点」、$Q_{ij}$ は「二つの選択肢を同時に選んだときの得点」と読めます。たとえば、ある道路を使う車が増えると混雑する、二つの商品を近くに置くと売れ方が変わる、同じ時間に同じ通路へAGVが入ると衝突する、という関係は $x_i{x}_j$ の項で表せます。

制約は、もっとも素朴には二乗ペナルティで入れます。

$$E(\mathit{x})=E_0(\mathit{x})+\sum _a{\lambda }_a{\{g_a(\mathit{x})\}}^2$$

$E_0$ は本来小さくしたい評価関数、$g_a(\mathit{x})=0$ は満たしたい制約、${\lambda }_a$ は制約違反の重みです。たとえば候補の中からちょうど一つ選ぶ制約は

$${(\sum _{\mu \in M_i}{x}_{\mu i}-1)}^2$$

と書けます。$i$ は車、商品、仕事、患者、移動体などの対象を表し、$\mu$ はその対象に用意された候補を表します。この「対象ごとに一つ選ぶ」という構造は、量子アニーリング応用のほとんどに現れます。

ただし、二乗ペナルティは無料ではありません。$\sum _{\mu }{x}_{\mu i}$ を二乗すると、同じ対象に属する候補どうしがすべて結合します。制約が大きくなるほど、QUBOは密になります。ここで大関法、半分ずらしの制約、列生成、Benders分解、クラスタ分解などが必要になります。詳しい分解法は 補遺: 制約つきQUBOの分解法 にまとめます。

2. 道路交通: 混雑を二次式で見る

量子アニーリングの実応用を語るとき、道路交通最適化はひとつの原点です。VolkswagenとD-Waveの共同研究では、北京の実タクシーデータをもとに、複数の車に対して複数の候補経路を用意し、どの車がどの経路を選ぶかをQUBOで決めました。車 $i$ が候補経路 $\mu$ を選ぶとき $x_{\mu i}=1$、選ばないとき $0$ とします。

道路リンクを $e$ とし、候補経路 $\mu$ がリンク $e$ を含むなら $C_{e\mu i}=1$、含まないなら $0$ とします。そのリンクを通る車の数は

$$C_e(\mathit{x})=\sum _i\sum _{\mu \in M_i}{C}_{e\mu i}{x}_{\mu i}$$

です。混雑を避けたいなら、同じリンクに車が集まるほど大きなコストを与えればよいです。最初の自然な書き方は

$$E_{\mathrm{t}\mathrm{r}\mathrm{a}\mathrm{f}\mathrm{f}\mathrm{i}\mathrm{c}}(\mathit{x})=\sum _e{\{C_e(\mathit{x})\}}^2+\frac{\lambda }{2}\sum _i{(\sum _{\mu \in M_i}{x}_{\mu i}-1)}^2$$

です。第一項は混雑、第二項は各車がちょうど一つの経路を選ぶ制約です。$C_e^2$ を展開すると、同じ道路を共有する二台の車の組に $x_{\mu i}{x}_{\nu j}$ の相互作用が出ます。つまり、混雑とは「同じ場所を同時に使う選択肢どうしの相互作用」なのです。

この式には、あとで大事な反省点が見つかります。$C_e^2$ は混雑の二体相互作用だけでなく、$C_e$ という一次の距離コストも含んでいます。なぜなら

$$C_e^2={(C_e-\frac{1}{2})}^2+C_e-\frac{1}{4}$$

だからです。右辺を見ると、$C_e^2$ は「リンクを誰かが使うだけで増える一次コスト」と「複数台が重なると増える混雑コスト」の和になっています。道路の長さや移動距離も同時に最小化したいなら一次コストが入ってよいですが、純粋に混雑を見たいなら、一次コストは余計です。

そこで、混雑そのものを取り出す自然な形は

$$E_{\mathrm{c}\mathrm{o}\mathrm{n}\mathrm{g}}(\mathit{x})=\sum _e{\{C_e(\mathit{x})-\frac{1}{2}\}}^2+\frac{\lambda }{2}\sum _i{(\sum _{\mu \in M_i}{x}_{\mu i}-1)}^2$$

です。$C_e=0$ と $C_e=1$ の値が同じになり、$C_e\ge 2$ になったときだけ追加コストが出ます。これは「道路は空でも一台でも構わないが、二台以上が重なると混雑です」という読み方に合っています。実応用のQUBOでは、この半分だけずらす発想が非常に効きます。衝突、容量、同時使用の禁止など、多くの制約がこの形で整理できます。

道路交通の研究が重要なのは、単に量子アニーラで交通を解いたからではありません。現場の目的関数をQUBOにしたとき、そこに何が余分に入ってしまうかまで見直したからです。QUBOは書ければ終わりではなく、書いた式を読む必要があります。

3. 津波避難と経路選択: 一人一経路から大規模避難へ

災害時の避難経路計画も、道路交通QUBOの自然な発展です。避難者、車両、またはグループ $i$ に複数の候補経路 $j\in R_i$ を用意し、$x_{ij}=1$ なら経路 $j$ を選ぶとします。移動時間や距離を $d_{ij}$、二つの経路が同じ道路を同じ時間に使うことで生じる干渉を $c_{ijkl}$ とすれば、基本形は

$$F(\mathit{x})=\alpha \sum _i\sum _{j\in R_i}{d}_{ij}{x}_{ij}+\sum _{i\ne k}\sum _{j\in R_i}\sum _{l\in R_k}{c}_{ijkl}{x}_{ij}{x}_{kl}$$

です。第一項は個々の移動負担、第二項は同時利用による混雑です。各対象が一つの経路を選ぶ制約は

$$E_{\mathrm{e}\mathrm{v}\mathrm{a}\mathrm{c}}(\mathit{x})=F(\mathit{x})+\beta \sum _i{(\sum _{j\in R_i}{x}_{ij}-1)}^2$$

となります。

避難計画で難しいのは、最短経路を一人ずつ選べばよいわけではない点です。全員が同じ橋、同じ交差点、同じ道路へ集まれば、個々には近道でも全体としては遅くなります。QUBOの二次項は、この「自分にとってよい道」と「みんなが選ぶと悪くなる道」の緊張関係を表します。

大規模避難では、全員の候補経路を一度にQUBO化すると変数数と結合数が膨らみます。そこで、地域や車両群をクラスタに分け、各クラスタ内でいくつかの候補解を作り、最後にクラスタごとの候補解を選ぶという二段階の方法が有効になります。クラスタ $m$ の候補解 $n$ を選ぶ変数を $s_{mn}$ とすれば、

$$H(\mathit{s})=\sum _{m,n}{A}_{mn}{s}_{mn}+\sum _{m\ne p}\sum _{n,q}{B}_{mnpq}{s}_{mn}{s}_{pq}+\gamma \sum _m{(\sum _n{s}_{mn}-1)}^2$$

のような上位QUBOになります。$A_{mn}$ はクラスタ内の候補解のコスト、$B_{mnpq}$ は二つのクラスタ解を同時に選んだときの干渉です。量子アニーラのサンプルが一つの最適解だけでなく多様な低エネルギー解を出すことは、この二段階化で役に立ちます。候補解の引き出しが増えるからです。

4. 二重制約と大関法: 全結合ペナルティを避ける

交通、避難、配送、割当では、「各対象が一つを選ぶ」と「各時刻や各場所に一つだけ入る」のように、二種類の制約が同時に出ます。典型的には、$x_{it}$ を「対象 $i$ が位置、時刻、またはスロット $t$ を選ぶ」変数として、

$$f(\mathit{x})=\sum _{i,t}{h}_{it}{x}_{it}+\frac{\lambda }{2}\sum _t{(\sum _i{x}_{it}-1)}^2+\frac{\lambda }{2}\sum _i{(\sum _t{x}_{it}-1)}^2$$

のような形になります。二つの一対一制約は分かりやすいのですが、展開すると同じ行、同じ列に属する変数がすべて結合します。変数が $N=L^2$ 個なら、埋め込みの負担は急に重くなります。

大関法は、この二乗ペナルティを補助場、つまりラグランジュ乗数の反復更新へ置き換えます。固定した補助場 $\mathit{\nu }$ のもとで

$$H_{\mathrm{e}\mathrm{f}\mathrm{f}}(\mathit{x};\mathit{\nu })=f_0(\mathit{x})-\sum _k{\nu }_k{F}_k(\mathit{x})$$

をサンプラで解き、得られたサンプルの平均が制約値に近づくように ${\nu }_k$ を更新します。二乗で全結合を作る代わりに、補助場で局所磁場を動かすわけです。ハバード・ストラトノビッチ変換から見た導出は 補遺: ハバード・ストラトノビッチ変換と大関法 にあります。不等式制約まで含めた実装上の整理は 補遺: 制約つきQUBOの分解法 で扱います。

この考え方は、量子アニーラの物理結合が限られているという弱点から生まれています。しかし、見方を変えると、制約を罰金で押しつぶすのではなく、双対変数を学習しながら問題を解く方法です。現場の大きな問題へ進むほど、この発想は重要になります。

5. 商品陳列: 売上と多様性のせめぎ合い

ECサイトの商品表示や棚割りでは、各商品をどの場所へ置くかを決めます。商品 $i$ を位置 $j$ へ置くなら $x_{ij}=1$ とします。人気や売上期待値を $s_{ij}$ とすれば、単純には

$$S(\mathit{x})=\sum _i\sum _j{s}_{ij}{x}_{ij}$$

を大きくしたくなります。しかし、人気順に並べるだけでは、似た商品ばかりが近くに集まったり、利用者が比較しにくい表示になったりします。そこで、商品どうしの類似度を $f_{i{i}^{\prime }}$、表示位置どうしの近さを $d_{j{j}^{\prime }}$ として、多様性の項を

$$D(\mathit{x})=-\sum _{i,i^{\prime }}\sum _{j,j^{\prime }}{f}_{i{i}^{\prime }}{d}_{j{j}^{\prime }}{x}_{ij}{x}_{i^{\prime }{j}^{\prime }}$$

のように入れます。似た商品を近い位置に置くほどコストを上げる、あるいは多様性の得点を下げる項です。商品は一つの場所にだけ置かれ、各場所には一つの商品だけ置くので、QUBOは

$$E_{\mathrm{l}\mathrm{i}\mathrm{s}\mathrm{t}}(\mathit{x})=-S(\mathit{x})-wD(\mathit{x})+M\sum _i{(\sum _j{x}_{ij}-1)}^2+M\sum _j{(\sum _i{x}_{ij}-1)}^2$$

と書けます。これは二次割当問題の一種です。売れるものを前に出すだけなら線形割当ですが、商品の組合せ効果や近接効果を入れると二次になります。

ここで量子アニーリングが担う役割は、ランキングを一つ返すだけではありません。低エネルギーの複数候補を返すことで、運用者が「売上寄り」「多様性寄り」「在庫消化寄り」などの候補を比較できます。現場では、最適解の一点よりも、説明できる良い候補の束が欲しい場面が多いのです。

6. 無人搬送車: 工場の通路を時間つきで予約する

DENSOの工場内AGV制御では、複数の無人搬送車が同じ工場内を走ります。各AGVには、現在位置と目的地からいくつかの候補経路が用意されます。AGV $i$ が候補経路 $\mu$ を選ぶなら $q_{\mu i}=1$ とします。候補経路の移動距離または進行度を $d_{\mu }$ とすれば、短い時間窓の中でできるだけ前へ進める目的は

$$f_0(\mathit{q})=-\sum _i\sum _{\mu \in M_i}{d}_{\mu }{q}_{\mu i}$$

と書けます。マイナス符号は、移動距離を大きくしたいからです。

衝突を避けるには、通路や交差点を表すエッジ $e$ と時刻 $t$ について、候補経路 $\mu$ がその時刻にそのエッジを使うかを $F_{\mu te}$ で表します。その場所・時刻に入るAGV数は

$$n_{e,t}(\mathit{q})=\sum _i\sum _{\mu \in M_i}{F}_{\mu te}{q}_{\mu i}$$

です。各AGVが一つの経路を選ぶ制約と、同じ場所・時刻に二台以上入らない制約を入れると、

$$f_{\mathrm{A}\mathrm{G}\mathrm{V}}(\mathit{q})=-\sum _i\sum _{\mu \in M_i}{d}_{\mu }{q}_{\mu i}+{\lambda }_1\sum _i{(\sum _{\mu \in M_i}{q}_{\mu i}-1)}^2+{\lambda }_2\sum _{e,t}{(n_{e,t}(\mathit{q})-\frac{1}{2})}^2$$

となります。最後の半分ずらしは、$n_{e,t}=0$ と $n_{e,t}=1$ を同じように許し、$n_{e,t}\ge 2$ を罰するためのものです。初期の定式化では制約を二乗で強く入れ、得られたサンプルから実行可能なものを選びました。リバースアニーリングを使う研究では、貪欲法などで作った実行可能解から出発し、量子揺らぎで近傍を探索します。

リバースアニーリング版では、すでにある候補解の周辺を調べるため、目的関数は残り距離 $d_{ij}^{\ast }$ を使って

$$f_{\mathrm{R}\mathrm{A}}(\mathit{q})=\sum _i\sum _j{d}_{ij}^{\ast }{q}_{ij}+a\sum _i{(\sum _j{q}_{ij}-1)}^2+b\sum _{e,t}{(\sum _i\sum _j{F}_{ijte}{q}_{ij}-\frac{1}{2})}^2$$

と書けます。初期解が悪くないとき、全空間を最初から探すよりも、その周辺を揺らして良い改善を探す方が自然です。これが工場内AGVとリバースアニーリングの相性のよさです。

大規模な商用AGVでは、1000台級の搬送車、動的に変わる優先度、通信遅延、リアルタイム性が入ります。ここでは、候補経路生成、優先度つきのQUBO、クラスタ分割、サンプラ、後処理が一体になります。量子アニーリングは、巨大な工場制御を単独で丸呑みするのではなく、短時間で多様な候補を出す部品として使われます。

7. 通信路復号: もっともありそうな誤りを探す

通信路復号では、受信語から送信語を推定します。線形符号なら、誤りベクトル $\mathit{e}$ とパリティ検査行列 $H$ に対して、シンドローム条件

$$H\mathit{e}=\mathit{s}(\mathrm{mod}2)$$

を満たす誤りの中から、通信路の確率に照らしてもっともありそうなものを選びます。ビット反転確率が小さいなら、誤りの重みを小さくしたいので、

$$E_{\mathrm{d}\mathrm{e}\mathrm{c}}(\mathit{e})=\sum _i{w}_i{e}_i+\lambda \sum _a{P}_a({\mathit{e}}_{\partial a},s_a)$$

のように書けます。$P_a$ は検査式 $a$ が破れたときに1、満たされたときに0になるペナルティです。高重みのパリティ検査はそのままでは高次項になりますが、LDPC符号のように疎な検査行列を使う、補助変数で二次化する、あるいは検査行列側の表現を工夫することでQUBOへ落とします。

この応用では、QUBO化のよさが符号表現に強く依存します。生成行列で書くか、パリティ検査行列で書くかによって、必要な変数数と結合数が変わります。実機へ載せるときには、数学的に同じ復号問題でも、どの表現を選ぶかが性能を決めます。

8. 鉄鋼製造スケジューリング: 二つの工程をつなぐ

鉄鋼製造では、複数の製品が複数の工程を通ります。前工程と後工程の順序、開始時刻、段取り替え、納期、設備容量が絡みます。製品 $i$ が工程 $p$ のスロット $s$ に入るなら $x_{p,i,s}=1$ とします。各製品は各工程で一つのスロットを取り、各スロットには一つの製品が入ります。

工程 $p$ における製品 $i$ の時刻を

$$t_{p,i}=\sum _{s}s{x}_{p,i,s}$$

とします。前工程が後工程より遅れてはいけないので、たとえば $t_{1,i}\le t_{2,i}$ が必要です。さらに、隣り合うスロットに置かれる製品の組合せによって段取り替えコストが変わります。

工程内の隣接コストを

$${\mathrm{C}\mathrm{o}\mathrm{s}\mathrm{t}}_{g,p}=\sum _{i,j}\sum _{s=1}^{N-1}{L}_{pij}{x}_{p,i,s}{x}_{p,j,s+1}$$

と書けば、順序づけの相互作用がQUBOの二次項になります。全体のQUBOには、工程ごとの割当制約、納期コスト、段取り替えコスト、工程間順序制約が入ります。

この問題は一気に解くと大きすぎるため、ラグランジュ分解協調法が使われます。工程間の結合制約をラグランジュ乗数 ${\lambda }_i$ で緩めると、

$$E_{\mathrm{L}\mathrm{D}\mathrm{C}}=\mathrm{C}\mathrm{o}\mathrm{s}\mathrm{t}+\rho \mathrm{P}\mathrm{e}\mathrm{n}\mathrm{a}\mathrm{l}\mathrm{t}\mathrm{y}+\sum _i{\lambda }_i(t_{1,i}-t_{2,i})$$

となり、工程1と工程2の問題へ分けられます。片方では $+{\lambda }_i{t}_{1,i}$、もう片方では $-{\lambda }_i{t}_{2,i}$ が局所的な重みとして働きます。解いた結果、前工程が後工程を追い越していれば ${\lambda }_i$ を更新します。制約を巨大な二乗ペナルティで入れるのではなく、二つの小さなQUBOを対話させるわけです。

9. ブラックボックス最適化と材料探索

学習から候補生成までの詳しい手順と研究史は、補遺: BOCSとFMによるブラックボックス最適化にまとめました。BOCSの事後分布による探索、FMの低ランク表現、SA・QAとの接続を計算例から追えます。

材料探索や分子設計では、候補を評価する実験やシミュレーションが高価です。目的関数の形が明示的に分からないので、直接QUBOを書くことができません。そこで、観測済みデータから二次の代理モデルを学習します。

$$\hat{f}(\mathit{x})=c+\sum _i{h}_i{x}_i+\sum _{i<j}{J}_{ij}{x}_i{x}_j$$

$\mathit{x}$ は材料の特徴、置換元素の有無、分子記述子、あるいは探索空間のビット表現です。代理モデル $\hat{f}$ をQUBOとして量子アニーラや古典サンプラに渡し、低い値を持つ候補を複数取り出します。その候補を本物の評価関数で測り、データに追加し、代理モデルを更新します。

これはブラックボックス最適化の一種です。量子アニーラが本質的に解いているのは代理モデルであり、実験室の物性値そのものではありません。しかし、巨大な候補空間から「次に試す価値のある候補」を出すサンプラーとして働きます。D-Wave実機とシミュレーテッドアニーリングを比較すると、単純に量子がいつも勝つわけではありません。それでも、低エネルギー候補の多様性、制約つき候補の生成、特徴量の組合せ探索という観点で、QUBO化は材料探索の設計言語になります。

仮想スクリーニングでは、化学空間が ${10}^{60}$ を超えるほど巨大になります。すべてを評価することはできないので、特徴量の重要性をサンプルから読み取り、探索空間を1%未満に絞るような使い方が考えられます。ここでも、量子アニーリングは「正解を一発で当てる装置」ではなく、広すぎる空間を読める大きさへ圧縮する道具として現れます。

10. 機械学習: データを選び、分布を学び、行列を分ける

量子アニーリング応用は、機械学習にも広がります。代表的な入り口は三つあります。第一に、学習データを選ぶ問題です。第二に、Boltzmann型モデルをサンプルする問題です。第三に、行列分解や特徴選択の二値変数を最適化する問題です。

誤ラベルを含む訓練データの選別では、訓練データ $i$ を採用するなら $x_i=1$、捨てるなら $0$ とします。採用集合で学習したモデルの検証損失を直接QUBOにすることは難しいので、低次の代理モデル

$$\hat{L}(\mathit{x})=c+\sum _i{a}_i{x}_i+\sum _{i<j}{b}_{ij}{x}_i{x}_j$$

を更新しながら、低損失を与えるデータ選択をサンプルします。誤ったラベルを持つデータは、良い採用集合の中に現れにくくなります。ここで大事なのは、最適な採用集合を一つ決めることだけではありません。複数の低損失サンプルから、どのデータが安定して選ばれるかを見ることです。

制限Boltzmann機械では、可視変数 $\mathit{v}$ と隠れ変数 $\mathit{h}$ に対して

$$E(\mathit{v},\mathit{h})=-\sum _{i,j}{w}_{ij}{v}_i{h}_j-\sum _i{b}_i{v}_i-\sum _j{c}_j{h}_j$$

というエネルギーを使います。学習にはモデル分布からのサンプルが必要です。量子アニーラをサンプラとして使うと、単なる最適化ではなく、低エネルギー状態の分布を得る問題になります。ただし実機の有効温度、スケール、局所磁場のずれを補正しないと、期待したBoltzmann分布にはなりません。オンラインキャリブレーションは、そのずれを推定しながら学習に使うための方法です。 フリーズ現象と有効温度の意味、Kimらの摂動論による温度・結合係数の予測は、ボルツマンマシンの補遺で導出します。

非負値二値行列分解では、非負行列 $V$ を

$$V\simeq WH$$

と分解します。$W$ は非負実数、$H$ は0-1行列です。目的関数は

$$\underset{W,H}{min}{\Vert V-WH\Vert }_{\mathrm{F}}^2,H\in \{0,1{\}}^{k\times n}$$

です。$W$ を固定すると、$H$ の各列 $\mathit{h}$ は独立に

$$E_j(\mathit{h})={\mathit{h}}^{\mathsf{T}}{W}^{\mathsf{T}}W\mathit{h}-2{\mathit{v}}_j^{\mathsf{T}}W\mathit{h}$$

を最小化するQUBOになります。顔画像、音源分離、特徴抽出では、二値の基底選択が解釈性を持ちます。リバースアニーリングを使えば、連続緩和で得た候補から出発し、二値解の近傍を探索できます。

カーネル学習や特徴選択でも、どの特徴を使うかを0-1変数で表し、性能を二次代理モデルやサンプル分布で評価する形が現れます。機械学習応用での量子アニーリングは、最適化器であり、候補集合の構造を可視化するサンプラーでもあります。

11. 交通信号制御: 信号の相を選ぶ

交通信号制御では、各交差点 $i$ でどの信号相 $p$ を選ぶかを決めます。$x_{ip}=1$ なら交差点 $i$ が相 $p$ を採用します。各交差点では一つの相だけを選ぶので、

$$\sum _p{x}_{ip}=1$$

です。交通流シミュレータから、ある相を選んだときの待ち時間や、隣接交差点との組合せで生じる渋滞コストを推定すれば、

$$E_{\mathrm{s}\mathrm{i}\mathrm{g}\mathrm{n}\mathrm{a}\mathrm{l}}(\mathit{x})=\sum _{i,p}{a}_{ip}{x}_{ip}+\sum _{(i,j)}\sum _{p,q}{b}_{ijpq}{x}_{ip}{x}_{jq}+M\sum _i{(\sum _p{x}_{ip}-1)}^2$$

と書けます。$a_{ip}$ は単独の待ち時間、$b_{ijpq}$ は隣接信号の相の組合せが作る流れの悪さです。

実都市の交通信号では、信号制御だけでなく、車両流入、歩行者、右左折、制御周期、シミュレータの近似が効きます。量子アニーラが返すQUBO解は、SUMOのような交通流シミュレータで評価されます。結果は、量子アニーラが常に古典ソルバを上回るという単純な話ではありません。しかし、信号制御を局所相の組合せ問題として切り出すことで、都市交通の一部をサンプラで実験できる形にしました。

12. ジョブショップと資源割当: まず資源を決め、次に時刻を決める

ジョブショップスケジューリングでは、どの資源を使うか、どの順番で作業するか、いつ始めるかを同時に決める必要があります。すべてを一つのQUBOに詰めると巨大になります。そこで、資源割当部分をQUBO、時刻決定部分を混合整数計画として分ける方法が使われます。

顧客または注文を $c$、タスクを $i$、機械を $j$、レイヤまたは候補を $l$ とします。タスク $i$ を機械 $j,l$ に割り当てるなら $x_{cijl}=1$、その機械を使うなら $y_{cjl}=1$ とします。基本のQUBOは

$$H=H_{\mathrm{f}\mathrm{i}\mathrm{l}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{g}}+{\lambda }_{\mathrm{d}\mathrm{e}\mathrm{a}\mathrm{d}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{e}}{H}_{\mathrm{d}\mathrm{e}\mathrm{a}\mathrm{d}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{e}}+{\lambda }_{\mathrm{o}\mathrm{n}\mathrm{e}-\mathrm{h}\mathrm{o}\mathrm{t}}{H}_{\mathrm{o}\mathrm{n}\mathrm{e}-\mathrm{h}\mathrm{o}\mathrm{t}}+{\lambda }_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{a}\mathrm{c}\mathrm{i}\mathrm{t}\mathrm{y}}{H}_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{a}\mathrm{c}\mathrm{i}\mathrm{t}\mathrm{y}}+{\lambda }_{xy}{H}_{xy}$$

です。各項は、使う資源数を減らす、納期や類似タスクのズレを抑える、一つのタスクを一つの資源に割り当てる、機械容量を守る、$x$ と $y$ の整合性を保つ、という役割を持ちます。たとえば、

$$H_{\mathrm{o}\mathrm{n}\mathrm{e}-\mathrm{h}\mathrm{o}\mathrm{t}}=\sum _{c,i}{(\sum _{j,l}{x}_{cijl}-1)}^2$$

であり、容量制約は

$$H_{\mathrm{c}\mathrm{a}\mathrm{p}\mathrm{a}\mathrm{c}\mathrm{i}\mathrm{t}\mathrm{y}}=\sum _{c,j,l}{(\sum _i{w}_{ci}{x}_{cijl}-\alpha B_j{y}_{cjl})}^2$$

のように書けます。$x_{cijl}\le y_{cjl}$ という含意制約は、半分ずらしに似た形

$$H_{xy}=\sum _{c,i,j,l}{(y_{cjl}-x_{cijl}-\frac{1}{2})}^2$$

で扱えます。$x=1,y=0$ を強く罰し、$x=0,y=0$ と $x=0,y=1$ と $x=1,y=1$ の扱いを目的関数側へ委ねる書き方です。

資源割当QUBOで候補を出したあと、具体的な開始時刻は混合整数計画で詰めます。量子アニーリングは、連続時刻まで含む全問題を一気に解くのではなく、組合せ部分の候補生成を担当します。現場問題では、この分担がとても自然です。

13. マイクロモビリティ配車: 待ち需要と充電を同時に見る

電動キックボードや小型移動体の配車では、使い終わった車両をどのステーションへ戻すか、あるいはどの充電設備へ運ぶかを決めます。車両 $v_i$ をステーション $s_j$ へ送るなら ${\sigma }_{v_i,s_j}=1$、充電器 $c_k$ へ送るなら ${\sigma }_{v_i,c_k}=1$ とします。各車両は一つの行き先を選ぶので、

$$H_{A0}=\sum _i{(1-\sum _j{\sigma }_{v_i,s_j}-\sum _k{\sigma }_{v_i,c_k})}^2$$

です。各充電器に同時に割り当てられる車両数を制御するなら、

$$H_{A1}=\sum _k{(1-\sum _i{\sigma }_{v_i,c_k})}^2$$

を入れます。移動時間のコストは

$$H_{B0}=\frac{1}{t_{\mathrm{a}\mathrm{v}\mathrm{g}}}\{\sum _{i,j}t(x_{v_i},x_{d_i},x_{s_j}){\sigma }_{v_i,s_j}+\sum _{i,k}t(x_{v_i},x_{d_i},x_{c_k}){\sigma }_{v_i,c_k}\}$$

です。$x_{v_i}$ は車両位置、$x_{d_i}$ は目的地、$x_{s_j}$ と $x_{c_k}$ はステーションと充電器の位置です。

さらに、ステーション $s_j$ で見込まれる需要を ${\tau }_j$ とし、そこへ戻す車両数を合わせる項を

$$H_{B1}=\sum _j{({\tau }_j-\sum _i{\sigma }_{v_i,s_j})}^2$$

とします。全体は

$$H=H_{A0}+H_{A1}+B_0{H}_{B0}+B_1{H}_{B1}$$

です。需要 ${\tau }_j$ は、過去の利用データからBayes的に推定できます。たとえば目的地候補 $c$ の確率を $p(x_c)$ とし、距離や移動時間に応じてステーションへの寄与を重みづければ、

$$P(s_j)=\int d{x}_{c}p(x_c)\frac{\exp [-t(s_j,c)/{\theta }_c]}{\sum _l\exp [-t(s_l,c)/{\theta }_c]}$$

のような形になります。現場での「需要がありそうな場所へ戻す」という直観が、確率分布を介してQUBOの係数へ入ります。

14. Benders分解: QUBOと連続最適化をつなぐ

混合整数二次計画では、二値変数 $\mathit{x}$ と連続変数 $\mathit{y}$ が同時に出ます。量子アニーラは二値変数を扱うのが得意ですが、連続変数はそのままでは扱いにくいです。Benders分解では、二値変数を主問題、連続変数を副問題へ分けます。

$$\underset{\mathit{x},\mathit{y}}{min}{\mathit{x}}^{\mathsf{T}}C\mathit{x}+{\mathit{h}}^{\mathsf{T}}\mathit{y}\mathrm{s}.\mathrm{t}.A\mathit{x}+G\mathit{y}\le \mathit{b}$$

とします。$\mathit{x}$ を固定すれば、$\mathit{y}$ については線形計画になります。その双対解から、$\mathit{x}$ が満たすべき切除平面を作ります。主問題は $\mathit{x}$ と補助変数 $t$ について

$$\underset{\mathit{x},t}{min}{\mathit{x}}^{\mathsf{T}}C\mathit{x}+t$$

を解き、得られた $\mathit{x}$ に対して副問題を解き、必要なら切除平面を追加します。$t$ やスラックを二値展開すると、この主問題をQUBOとして扱えます。Benders分解の詳しい流れは 補遺: 制約つきQUBOの分解法 にまとめます。

この方法の良いところは、量子アニーラに全問題を押し付けない点です。二値の組合せだけを量子アニーラに渡し、連続最適化や切除平面は古典計算で厳密に扱います。現場の問題は、二値、整数、連続、論理制約が混ざります。ハイブリッド分解は、その混ざり方を尊重する方法です。

15. 年表として見る応用の広がり

年	話題	QUBOで見た主変数	解法上の要点
2017	道路交通最適化	車ごとの候補経路	混雑をリンク共有の二次項で表す
2019	商品陳列	商品と表示位置の割当	売上と多様性の二次割当
2019	DENSO AGV	AGVごとの候補経路	時間つき通路占有と衝突回避
2020	通信路復号	誤りビット	パリティ検査制約のQUBO化
2020	二重制約の大関法	割当変数と補助場	全結合ペナルティを避ける
2021	ブラックボックス最適化	材料・設計特徴ビット	二次代理モデルをサンプルする
2022	鉄鋼スケジューリング	工程・製品・スロット	ラグランジュ分解協調
2022	リバースアニーリングAGV	実行可能解近傍の経路選択	初期解から局所探索する
2023	仮想スクリーニング	化学特徴ビット	サンプルで探索空間を絞る
2023	列生成	制約つき二次0-1問題の候補列	価格付け問題をQUBOにする
2025	交通信号制御	交差点ごとの信号相	交通流シミュレータと組み合わせる
2025	誤ラベルデータ選別	訓練データ採否	検証損失の二次代理モデル
2025	RBM学習	可視・隠れビット	実機サンプルの有効温度を校正する
2025	NBMFとリバースアニーリング	行列分解の二値係数	連続緩和解から二値近傍探索へ
2025	商用AGV大規模制御	1000台級の経路候補	優先度、クラスタ分割、短時間サンプリング
2026	Benders分解	二値主変数	連続副問題から切除平面を作る
2026	ジョブショップ	資源割当と時刻	QUBOと混合整数計画を分担する
2026	音源分離	NBMFの二値係数	行列分解をQUBO列問題に分ける
2026	マイクロモビリティ配車	車両の行き先	需要推定、充電、待ち時間を同時に扱う
2026	避難経路計画	避難主体ごとの候補経路	多様なサンプルとクラスタ合成

この年表で見ると、応用は単純な「QUBOへ直して解く」から、より現実的な「QUBOをハイブリッドアルゴリズムの一部として使う」へ移っています。初期には、交通やAGVのように、候補の組合せを直接QUBOにしました。その後、大規模化、制約処理、サンプルの多様性、連続変数との分担、機械学習との接続が前面に出ています。

量子アニーリング応用の歴史は、量子アニーラが強くなった歴史であり、QUBOを書く側が賢くなった歴史でもあります。現場の問題は、きれいな二次式として待っていてくれません。だからこそ、定式化の工夫が研究になります。

16. QUBO設計で見るべきもの

応用論文を読むときは、次の四つを見ます。

第一に、何が0-1変数かです。車と経路、商品と位置、AGVと経路、データと採否、特徴量と採否、資源とタスク。変数の選び方が問題の骨格を決めます。

第二に、二次項が何を表すかです。混雑、衝突、類似商品の近接、隣接工程の段取り、信号の相互干渉、特徴量どうしの相乗効果。二次項は、現場の相互作用です。

第三に、制約をどう処理しているかです。普通の一対一制約なら二乗ペナルティでよい場合があります。容量制約、不等式、二重制約、連続変数が入るなら、大関法、半分ずらし、列生成、Benders分解、ラグランジュ分解、クラスタ分解が必要になります。

第四に、量子アニーラがどこで効いているかです。最適解を一つ見つけるのか、低エネルギー解をたくさん出すのか、既存解の近傍をリバースアニーリングで探すのか、代理モデルの候補を出すのでしょうか。ここを混同すると、実用性の議論がぼやけます。

量子アニーリングの応用でいちばん大事なのは、量子という言葉に問題を合わせることではありません。問題のほうをよく見て、どこに二値の選択があり、どこに相互作用があり、どこに譲れない制約があるかを見つけることです。その翻訳ができると、量子アニーラは現場の中で役割を持ち始めます。

参考文献

1. Florian Neukart, Gabriele Compostella, Christian Seidel, David von Dollen, Sheir Yarkoni, Bob Parney, “Traffic Flow Optimization Using a Quantum Annealer,” Frontiers in ICT 4, 29 (2017). 出典

2. Masayuki Ohzeki, “Reconsideration of optimization for reduction of traffic congestion,” arXiv:2406.05448 (2024). 出典

3. Masayuki Ohzeki, “Breaking limitation of quantum annealer in solving optimization problems under constraints,” Scientific Reports 10, 3126 (2020). 出典

4. Naoki Nishimura, Kotaro Tanahashi, Koji Suganuma, Masamichi J. Miyama, Masayuki Ohzeki, “Item Listing Optimization for E-Commerce Websites Based on Diversity,” Frontiers in Computer Science 1, 2 (2019). 出典

5. Masayuki Ohzeki, Akira Miki, Masamichi J. Miyama, Masayoshi Terabe, “Control of Automated Guided Vehicles Without Collision by Quantum Annealer and Digital Devices,” Frontiers in Computer Science 1, 9 (2019). 出典

6. Renichiro Haba, Masayuki Ohzeki, Kazuyuki Tanaka, “Travel time optimization on multi-AGV routing by reverse annealing,” arXiv:2204.11789 (2022). 出典

7. Kouki Yonaga, Masamichi Miyama, Masayuki Ohzeki, Koji Hirano, Hirokazu Kobayashi, Tetsuaki Kurokawa, “Quantum Optimization with Lagrangian Decomposition for Multiple-process Scheduling in Steel Manufacturing,” ISIJ International 62, 1874–1880 (2022). 出典

8. Ami S. Koshikawa, Masayuki Ohzeki, Tadashi Kadowaki, Kazuyuki Tanaka, “Benchmark Test of Black-box Optimization Using D-Wave Quantum Annealer,” Journal of the Physical Society of Japan 90, 064001 (2021). 出典

9. Mikiya Doi, Yoshihiro Nakao, Takuro Tanaka, Masami Sako, Masayuki Ohzeki, “Exploration of new chemical materials using black-box optimization with the D-wave quantum annealer,” Frontiers in Computer Science 5, 1286226 (2023). 出典

10. Takuro Tanaka, Masami Sako, Mahito Chiba, Chul Lee, Hyukgeun Cha, Masayuki Ohzeki, “Virtual Screening of Chemical Space Based on Quantum Annealing,” Journal of the Physical Society of Japan 92, 023001 (2023). 出典

11. Makoto Otsuka, Kento Kodama, Keisuke Morita, Masayuki Ohzeki, “Filtering out mislabeled training instances using black-box optimization and quantum annealing,” Scientific Reports 15, 37892 (2025). 出典

12. Renichiro Haba, Masayuki Ohzeki, Kazuyuki Tanaka, “Relaxation-assisted reverse annealing on nonnegative/binary matrix factorization,” PLOS One 20, e0323232 (2025). 出典

13. Reo Shikanai, Masayuki Ohzeki, Kazuyuki Tanaka, “Quadratic Unconstrained Binary Formulation for Traffic Signal Optimization on Real-World Maps,” Journal of the Physical Society of Japan 94, 024001 (2025). 出典

14. Takeru Goto, Masayuki Ohzeki, “Online Calibration Scheme for Training Restricted Boltzmann Machines with Quantum Annealing,” Journal of the Physical Society of Japan 94, 034002 (2025). 出典

15. Thinh Nguyen Quang, Kosuke Matsuyama, Keisuke Shimizu, Hiroki Sugano, Eiji Kurimoto, Masahiro Miki, et al., “Quantum annealing-based route optimization for commercial AGV operating systems in large-scale logistics warehouses,” Scientific Reports 15, 44047 (2025). 出典

16. Sota Hirama, Masayuki Ohzeki, “Efficient Algorithm for Binary Quadratic Problem by Column Generation and Quantum Annealing,” Journal of the Physical Society of Japan 92, 113002 (2023). 出典

17. Taisei Takabayashi, Naoki Maruyama, Takuma Yoshihara, Renichiro Haba, Masayuki Ohzeki, “Efficient Construction of Feasible Solutions in Column Generation Using Quantum Annealing,” Journal of the Physical Society of Japan 94, 094005 (2025). 出典

18. Takuma Yoshihara, Masayuki Ohzeki, “Accelerating Extended Benders Decomposition with Quantum-Classical Hybrid Solver,” Journal of the Physical Society of Japan 95, 033001 (2026). 出典

19. Kenta Sawamura, Kensuke Araki, Naoki Maruyama, Renichiro Haba, Masayuki Ohzeki, “Quantum-Classical Hybrid Algorithm Using Quantum Annealing for Multi-Objective Job Shop Scheduling,” Journal of the Physical Society of Japan 95, 054002 (2026). 出典

20. Takuto Sasajima, Ami S. Koshikawa, Masayuki Ohzeki, “Audio Source Separation via Nonnegative Binary Matrix Factorization Using D-Wave Quantum Annealer,” Journal of the Physical Society of Japan 95, 064006 (2026). 出典

21. Yasushi Hasegawa, Hiroki Oshiyama, Masayuki Ohzeki, “Kernel Learning by Quantum Annealer,” Journal of the Physical Society of Japan 95, 074002 (2026). 出典

22. Takeru Goto, Masayuki Ohzeki, “Micro-mobility dispatch optimization via quantum annealing incorporating historical data,” Scientific Reports 16, 24054 (2026). 出典

23. Reo Shikanai, Renichiro Haba, Yusuke Okazaki, Kazumichi Matsumoto, Masayuki Ohzeki, “Large-scale evacuation route optimization leveraging sampling diversity in quantum annealing,” Scientific Reports 16, 22108 (2026). 出典