はじめに

AIREV株式会社の上里です。
量子アニーリングはQUBO（二次制約なし二値最適化）の形で定式化された組み合わせ最適化問題を量子コンピューティングによって解く手法です。したがって、目的の問題をQUBOの形に落とし込むことができれば量子アニーリングを適用することができます。D-waveの量子アニーリングのクラウドサービスLeapで数独を解いてみるでもQUBOと同等なBinary Quadratic Model（BQM）の形で問題の定式化を行いました。
しかし、QUBOやBQMでは直接的に制約を表現できない（目的関数にペナルティ項を設ける形になる）、変数のとる値が2値（バイナリ）でなければならない、などの特徴から問題の性質によっては直感と異なる形での定式化が必要になることがあります。
これに対して、Leapでは問題の性質に合わせてより直感に近い形で定式化を行えるよう、複数のモデルによる定式化をサポートしています。今回はその中から直接制約を表現することが可能なConstrained Quadratic Model（CQM）を紹介し、実際にCQMで定式化した問題をLeapを利用して解いてみます。

QUBOやBQMでの制約の表現方法

組み合わせ最適化問題は、与えられた制約のもとで目的の指標を最も良くするような変数の値の組み合わせを求める問題です。例えば容量の上限を超えずに最も価値の高い品物の組み合わせを選ぶナップサック問題では「容量の上限を超えず」の部分が制約にあたります。

QUBOやBQMで問題の定式化を行う場合は制約を直接表現することができないため、目的関数内に制約としての働きをする項（ペナルティ項）を設ける形で制約を表現します。QUBOで数独の問題を定式化する例はこちらの記事でご参照いただけます。制約の表現方法が少々難解に感じる方もいらっしゃるのではないでしょうか。

Constrained Quadratic Model（CQM）

CQMでは問題の制約を直接表現することができます。具体的には、 $C_{eq.}$ 個の等式制約と $C_{ineq.}$ 個の不等式制約、そして最小化する目的関数によって問題を定式化します。設定した制約のもとで目的関数を最小化する $N$ 個の変数 $\{x_i\}_{i=1,...,N}$ の値の組み合わせが問題の解となります。変数は $\{0, 1\}$ の2値を取るものと整数値を取るものを定義することができます。

目的関数：

\begin{aligned} \sum_{i} a_i x_i + \sum_{i<j} b_{ij} x_i x_j + c, \end{aligned}

制約：

\begin{aligned} \sum_{i} a_i^{(d)} x_i + \sum_{i<j} b_{ij}^{(d)} x_i x_j + c^{(d)} \leq 0, \qquad d = 1,...,C_{ineq.}, \end{aligned}

\begin{aligned} \sum_{i} a_i^{(e)} x_i + \sum_{i<j} b_{ij}^{(e)} x_i x_j + c^{(e)} = 0, \qquad e = 1,...,C_{eq.}, \end{aligned}

（ $a_i, b_{ij}, c$ は実数です。）

Leapを利用するためのSDKであるOcean SDKでは非常に直感的なインターフェースが提供されているため、上記の式を詳細に理解していなくてもSDKを利用した定式化を行うことは可能かと思います。ここで把握しておくとよい点は下記の2点です。

$\{0, 1\}$ の2値をとる変数と整数値をとる変数を設定できる
不等式制約・等式制約の両方がサポートされている

LeapでCQMを試してみる

早速Leapを利用してCQMで定式化した問題を解いてみましょう。LeapやOcean SDKの初期設定についてはこの記事の対象外となりますので、必要に応じて公式ドキュメントや解説記事等を参照してください。Ocean SDKはPythonパッケージとして配布されていますので、以下ではPythonでコーディングしていきます。

環境

この記事のコードを動作させるためにはOcean SDKが必要になります。インストール方法は公式ドキュメントをご参照ください。
CQMはOcean SDKのバージョン4.0.0以降でサポートされています。
この記事のコードはPython 3.9.1、Ocean SDK 4.2.0での動作を確認済みです。

問題設定

今回はビンパッキング問題を例題としたいと思います。ビンパッキング問題は、与えられた複数のアイテムをビンに詰めるとき、どのアイテムの組み合わせで詰めると最小のビンの本数で全てのアイテムを詰めることができるかを求める問題です。

今回は4つのアイテムがそれぞれ $0.9, 0.1, 0.7, 0.2$ の重みを持っており、1本のビンに $1.0$ までの重みのアイテムを詰められることとします。

アイテムの重みのリストを weights 、ビンのキャパシティを capacity とします。

weights = [.9, .1, .7, .2]
capacity = 1

変数

定式化をしていきます。まずは必要な変数の用意です。変数のクラスはOcean SDKの dimod モジュールで提供されています。

import dimod

# ビンの利用状況。それぞれの変数が各ビンの利用/不利用を表す。
y = [dimod.Binary(f'y_{j}') for j in range(len(weights))]
# それぞれの変数がどのアイテムがどのビンに入っているかを示す。jがビンのID、iがアイテムのID。
x = [[dimod.Binary(f'x_{i}_{j}') for j in range(len(weights))]
　　for i in range(len(weights))]

y はそれぞれのビンを利用しているか否かを示す2値の変数を要素としてもつリストです。十分なビンの数を用意したいので、アイテム数分のビン（変数）を用意しておきます。たとえば0番目のビンと2番目のビンが利用されている場合、 y は [1, 0, 1, 0] となります。
x は i 番目アイテムが j 番目のビンに入っているか否かを示す2値の変数を要素としてもつ2次元のリストです。i 番目のアイテムが j 番目のビンに詰められている場合は x[i][j] の値が 1 になります。該当しないアイテムとビンの組み合わせでは値が 0 になります。

今回の例では2値の変数のみを利用していますが、整数の変数は dimod.Integer クラスのインスタンスとして表現することができます。整数の変数の利用例はD-waveのサンプルをご参照ください。

目的関数

次は目的関数です。利用するビンの数を最小化したいので、下記のようになります。

cqm = dimod.ConstrainedQuadraticModel()
cqm.set_objective(sum(y))

CQMのインスタンスを dimod モジュールから取得し、Pythonの組み込み関数 sum を利用して y の要素の合計を最小化するようにCQMに目的関数をセットします。

制約

次に制約を設定していきます。まず、それぞれのアイテムが必ず1本のビンにだけ詰められるように制約を設定します。制約を設定する際はCQMインスタンスの add_constraint メソッドを利用できます。

for i in range(len(weights)):
    cqm.add_constraint(sum(x[i]) == 1, label=f'item_placing_{i}')

それぞれのアイテムに対して x[i] の要素の合計が 1 になるような等式制約を設定しました。
さらにそれぞれのビンのキャパシティを超えないように不等式制約を追加します。

for j in range(len(weights)):
    cqm.add_constraint(
        sum(weights[i] * x[i][j] for i in range(len(weights))) - y[j] * capacity <= 0,
        label=f'capacity_bin_{j}')

sum(weights[i] * x[i][j] for i in range(len(weights))) では j 番目のビンに詰められているアイテムの重みの合計を算出しています。
- y[j] * capacity としている部分の意図が少々分かりづらいかと思います。この部分は j 番目のビンに1つでもアイテムが詰められるときには必ず y[j] が 1 となるような工夫になっています。 j 番目のビンにアイテムが1つでも詰められていれば第一項の sum(...) が正の値となります。このとき y[j] が 0 であると、左辺全体の sum(...) - y[j] * capacity の値も正の値となってしまい、制約を満たすことができないことがわかると思います。

以上で定式化は完了です。

ハイブリッドソルバで解いてみる

CQMで定式化された問題を解く計算機として、LeapではCQM専用のハイブリッドソルバが用意されています。ハイブリッドソルバとは、古典コンピュータのCPU・GPUでの演算と量子アニーリングを併用して組み合わせ最適化問題を解くコンピュータです。今回はCQM専用のハイブリッドソルバを利用して、上述のビンパッキング問題を解いてみます。

from dwave.system import LeapHybridCQMSampler

sampler = LeapHybridCQMSampler()
sampleset = sampler.sample_cqm(cqm, label='CQM example')

# 制約を満たす解を取得
feasible_sample = sampleset.filter(lambda sample: sample.is_feasible).first

print(feasible_sample.energy) # 2.0
print(feasible_sample.sample)
# {'x_0_0': 0.0, 'x_0_1': 1.0, 'x_0_2': 0.0, 'x_0_3': 0.0,
#  'x_1_0': 0.0, 'x_1_1': 1.0, 'x_1_2': 0.0, 'x_1_3': 0.0,
#  'x_2_0': 0.0, 'x_2_1': 0.0, 'x_2_2': 0.0, 'x_2_3': 1.0,
#  'x_3_0': 0.0, 'x_3_1': 0.0, 'x_3_2': 0.0, 'x_3_3': 1.0,
#  'y_0': 0.0, 'y_1': 1.0, 'y_2': 0.0, 'y_3': 1.0}

dwave.system モジュールの LeapHybridCQMSampler は、Leapのハイブリッドソルバから得られた解のパターンを取得してくれます。解のパターンの中から設定した制約を満たしている解を取得し、得られた解の目的関数の値を energy プロパティから、また変数の値を sample プロパティから取得できます。
上記の結果では、0番目と1番目のアイテムを1番目のビンに、2番目と3番目のアイテムを3番目のビンに詰めることで、2本のビンを使って全てのアイテムを詰める解が求められました。

まとめ

今回は制約を直接的に設定することができる定式化モデルであるCQMを紹介し、ビンパッキング問題を例題として実際にLeapのハイブリッドソルバでCQMで定式化された問題の解を求めてみました。
離散値をとる変数を直感的に扱えるDiscrete Quadratic Model（DQM）についても後日別の記事で紹介しようと思います。

参考

（2021年11月26日アクセス）

https://github.com/dwave-examples
（2021年11月26日アクセス）

Discussion