ナップザック問題

重さと価値を持つアイテムの集合と、重量制限のあるナップサックが与えられたとき、ナップサック問題は、総重量が容量以内に収まるようにしつつ、総価値を最大化するアイテムの部分集合を選択することを目的とします。 \(w_i\)\(v_i\) (\(0 \leq i \leq n-1\)) をそれぞれアイテム \(i\) の重さと価値とします。 \(S \subset \{0, \cdots, n-1\}\) を選択されたアイテムの集合とします。

\[\begin{split} \mathrm{max} &\sum_{i \in S}v_i,\\ \mathrm{subject~to} &\sum_{i \in S}w_i \leq W. \end{split}\]

ここで \(W\) はナップサックの重量容量です。

QUBO定式化

この問題をQUBOとして定式化するために、\(n\) 個のバイナリ変数 \(x_i \in \{0,~1\}\) (\(0 \leq i \leq n-1\)) の集合 \(X\) を導入します。ここで、アイテム \(i\) が選択されるのは \(x_i=1\) のときかつそのときに限ります。

上記の定式化は次のように書き換えられます:

\[\begin{split} \mathrm{max} &\sum_{i=0}^{n-1}v_ix_i,\\ \mathrm{subject~to} &\sum_{i=0}^{n-1}w_ix_i \leq W. \end{split}\]

QUBO++プログラム

制約はHi-QUBOが提供する 範囲演算子 を用いて表現できます。 結果として得られるQUBO目的関数は次のように定義されます:

PyQBPPプログラム

制約はpyQBPPが提供する 範囲演算子 を用いて表現できます。 結果として得られるQUBO目的関数は次のように定義されます:

\[ f(X) = -\sum_{i=0}^{n-1}v_ix_i + P \left( 0 \leq \sum_{i=0}^{n-1}w_ix_i \leq W \right) \]

QUBOソルバーは目的関数を最小化するため、元の最大化目的は符号を反転しています。 定数 \(P\) は制約を強制するための十分大きなペナルティパラメータです。

以下のプログラムは、Exhaustive Solverを用いて10個のアイテムのナップサック問題を解きます:

cutting-stock-program.cpp
#include <qbpp/qbpp.hpp>
#include <qbpp/easy_solver.hpp>

int main() {
  const int L = 60;
  const auto l = qbpp::array({13, 23, 8, 11});
  const auto c = qbpp::array({10, 4, 8, 6});
  const size_t N = l.size();
  const size_t M = 6;

  auto x = qbpp::var_int("x", M, N) == 0;
  for (size_t i = 0; i < M; i++) {
    for (size_t j = 0; j < N; j++) {
      x[i][j] = 0 <= qbpp::var_int() <= c[j];
    }
  }

  auto order_fulfilled_count = qbpp::vector_sum(x, 0);
  auto order_constraint = order_fulfilled_count == c;

  auto bar_length_used = qbpp::expr(M);
  for (size_t i = 0; i < M; i++) {
    bar_length_used[i] = qbpp::sum(x(i) * l);
  }
  auto bar_constraint = 0 <= bar_length_used <= L;

  auto f = qbpp::sum(order_constraint) + qbpp::sum(bar_constraint);
  f.simplify_as_binary();

  auto solver = qbpp::EasySolver(f);
  auto sol = solver.search({{"time_limit", 10.0}, {"target_energy", 0}});
  for (size_t i = 0; i < M; i++) {
    std::cout << "Bar " << i << ":  ";
    for (size_t j = 0; j < N; j++) {
      std::cout << sol(x[i][j]) << "  ";
    }
    std::cout << " used = " << sol(bar_length_used[i])
              << ", waste = " << L - sol(bar_length_used[i]) << std::endl;
  }
  for (size_t j = 0; j < N; j++) {
    std::cout << "Order " << j
              << " fulfilled = " << sol(order_fulfilled_count[j])
              << ", required = " << c[j] << std::endl;
  }
}
cutting-stock-program.py
import pyqbpp as qbpp

L = 60
l = qbpp.array([13, 23, 8, 11])
c = qbpp.array([10, 4, 8, 6])
N = len(l)
M = 6

# 棒 i から切り出される注文 j のピース数を表す整数変数 x[i][j] を生成
x = [[qbpp.var(between=(0, c[j])) for j in range(N)] for i in range(M)]

# 注文制約:各注文について合計ピース数が c[j] と等しくなければならない
order_fulfilled_count = []
order_constraint = 0
for j in range(N):
    col_sum = 0
    for i in range(M):
        col_sum += x[i][j]
    order_fulfilled_count.append(col_sum)
    order_constraint += (col_sum == c[j])

# 棒制約:各棒で使用される合計長は L を超えてはならない
bar_length_used = []
bar_constraint = 0
for i in range(M):
    used = 0
    for j in range(N):
        used += x[i][j] * l[j]
    bar_length_used.append(used)
    bar_constraint += (0 <= used) & (qbpp.same <= L)

f = order_constraint + bar_constraint
f.simplify_as_binary()

solver = qbpp.EasySolver(f)
sol = solver.search(time_limit=10.0, target_energy=0)

for i in range(M):
    pieces = "  ".join(str(sol(x[i][j])) for j in range(N))
    used = sol(bar_length_used[i])
    print(f"Bar {i}:  {pieces}   used = {used}, waste = {L - used}")

for j in range(N):
    fulfilled = sol(order_fulfilled_count[j])
    print(f"Order {j} fulfilled = {fulfilled}, required = {c[j]}")

このプログラムでは、式 constraintobjective を別々に構築し、ペナルティ係数 1000 を用いて最終的なQUBO式 f に結合しています。 次に、Exhaustive Solverf に適用し、すべての最適解を列挙します。

以下の出力は、エネルギー、制約値、目的関数値を含む最適解を示しています:

[Solution 0]
Energy = -480
Constraint  = 50
Objective  = 480
Item 3: weight = 5, value =  60
Item 5: weight = 15, value =  150
Item 6: weight = 12, value =  110
Item 9: weight = 18, value =  160
[Solution 1]
Energy = -480
Constraint  = 50
Objective  = 480
Item 3: weight = 5, value =  60
Item 4: weight = 8, value =  80
Item 6: weight = 12, value =  110
Item 7: weight = 7, value =  70
Item 9: weight = 18, value =  160

このインスタンスには2つの最適解があり、いずれも総価値 480 を達成しつつ、容量制約をちょうど満たしていることがわかります。