変数の固定と番兵

変数の配列の一部の値が、問題のデータとして初めから決まっていることがあります。 QUBO++ でこうした変数を固定する方法は 2 通りあります:

• 方法 1: 変数の配列のまま式を作り、qbpp::replace で固定する値に置き換える
• 方法 2: 変数の配列を式の配列に移し、固定する要素に定数を代入してから式を作る

このページでは、まず小さなパズルで 2 つの方法を比べます。 そのあと、配列の外側に固定値（番兵）を置く場合を、同じ 2 つの方法で比べます。

行和・列和が決まった 0/1 行列

$4 \times 4$ の 0/1 行列 $x$ の各行・各列の 1 の個数が決まっています:

\[\begin{aligned} \sum_{j=0}^{3} x_{i,j} &= R_i \quad (i = 0, 1, 2, 3), & R &= (2, 3, 2, 2)\\ \sum_{i=0}^{3} x_{i,j} &= C_j \quad (j = 0, 1, 2, 3), & C &= (2, 3, 2, 2) \end{aligned}\]

この条件だけでは行列は 51 通りありますが、3 つのマスの値 $x_{0,1} = 0$、$x_{1,0} = 0$、$x_{3,3} = 1$ が分かっていると、1 通りに決まります:

\[x = \begin{pmatrix} 1 & 0 & 1 & 0\\ 0 & 1 & 1 & 1\\ 1 & 1 & 0 & 0\\ 0 & 1 & 0 & 1 \end{pmatrix}\]

分かっている 3 つのマスを固定して、残りのマスを求めます。

方法 1: replace で固定する

次のプログラムは、変数の配列 x のままで行和・列和の制約 f を作り、qbpp::replace で分かっているマスを固定します。 分かっているマスが違う 2 つの問題を、同じ f から解いています:

#include <qbpp/qbpp.hpp>
#include <qbpp/exhaustive_solver.hpp>

int main() {
  auto R = qbpp::array({2, 3, 2, 2});
  auto C = qbpp::array({2, 3, 2, 2});
  auto x = qbpp::var("x", 4, 4);
  auto f = qbpp::cons(qbpp::vector_sum(x, 1) == R) +
           qbpp::cons(qbpp::vector_sum(x, 0) == C);
  std::cout << "f =" << std::endl << f << std::endl;

  qbpp::MapList ml1({{x[0][1], 0}, {x[1][0], 0}, {x[3][3], 1}});
  qbpp::MapList ml2({{x[0][0], 0}, {x[0][1], 0}, {x[3][3], 1}});
  for (const auto& ml : {ml1, ml2}) {
    auto g = qbpp::replace(f, ml);
    g.simplify_as_binary();
    std::cout << "g =" << std::endl << g << std::endl;
    auto sol = qbpp::ExhaustiveSolver(g).search();
    auto full_sol = qbpp::Sol(f).set(sol, ml);
    std::cout << "x = " << full_sol(x) << std::endl;
  }
}

このプログラムの出力は次のとおりです:

f =
cons(x[0][0] +x[0][1] +x[0][2] +x[0][3] == 2)
+cons(x[1][0] +x[1][1] +x[1][2] +x[1][3] == 3)
+cons(x[2][0] +x[2][1] +x[2][2] +x[2][3] == 2)
+cons(x[3][0] +x[3][1] +x[3][2] +x[3][3] == 2)
+cons(x[0][0] +x[1][0] +x[2][0] +x[3][0] == 2)
+cons(x[0][1] +x[1][1] +x[2][1] +x[3][1] == 3)
+cons(x[0][2] +x[1][2] +x[2][2] +x[3][2] == 2)
+cons(x[0][3] +x[1][3] +x[2][3] +x[3][3] == 2)
g =
cons(x[0][0] +x[0][2] +x[0][3] == 2)
+cons(x[1][1] +x[1][2] +x[1][3] == 3)
+cons(x[2][0] +x[2][1] +x[2][2] +x[2][3] == 2)
+cons(1 +x[3][0] +x[3][1] +x[3][2] == 2)
+cons(x[0][0] +x[2][0] +x[3][0] == 2)
+cons(x[1][1] +x[2][1] +x[3][1] == 3)
+cons(x[0][2] +x[1][2] +x[2][2] +x[3][2] == 2)
+cons(1 +x[0][3] +x[1][3] +x[2][3] == 2)
x = {{1,0,1,0},{0,1,1,1},{1,1,0,0},{0,1,0,1}}
g =
cons(x[0][2] +x[0][3] == 2)
+cons(x[1][0] +x[1][1] +x[1][2] +x[1][3] == 3)
+cons(x[2][0] +x[2][1] +x[2][2] +x[2][3] == 2)
+cons(1 +x[3][0] +x[3][1] +x[3][2] == 2)
+cons(x[1][0] +x[2][0] +x[3][0] == 2)
+cons(x[1][1] +x[2][1] +x[3][1] == 3)
+cons(x[0][2] +x[1][2] +x[2][2] +x[3][2] == 2)
+cons(1 +x[0][3] +x[1][3] +x[2][3] == 2)
x = {{0,0,1,1},{1,1,1,0},{1,1,0,0},{0,1,0,1}}

f には 16 個の変数がすべて現れます。 qbpp::replace で $0$ に固定した変数は項ごと消え、$1$ に固定した x[3][3] は定数 1 になっています。

g には固定した変数が現れないので、ソルバーの解 sol にも含まれません。 qbpp::Sol(f).set(sol, ml) で、元の式 f のすべての変数を持つ解に戻してから、行列 x 全体を読みます（置換関数）。

2 つ目の問題は、f を作り直さずに、対応表 ml2 で置き換えるだけで作れます。

方法 2: 式の配列に定数を代入して作る

方法 2 では、変数の配列 x を式の配列 y に移し、分かっているマスに定数を代入してから、y で制約を作ります:

#include <qbpp/qbpp.hpp>
#include <qbpp/exhaustive_solver.hpp>

int main() {
  auto R = qbpp::array({2, 3, 2, 2});
  auto C = qbpp::array({2, 3, 2, 2});
  auto x = qbpp::var("x", 4, 4);
  auto y = qbpp::expr(4, 4);
  y += x;
  y[0][1] = 0;
  y[1][0] = 0;
  y[3][3] = 1;
  std::cout << "y = " << y << std::endl;

  auto f = qbpp::cons(qbpp::vector_sum(y, 1) == R) +
           qbpp::cons(qbpp::vector_sum(y, 0) == C);
  f.simplify_as_binary();
  std::cout << "f =" << std::endl << f << std::endl;

  auto sol = qbpp::ExhaustiveSolver(f).search();
  std::cout << "x = " << sol(y) << std::endl;
}

qbpp::expr(4, 4) はすべての要素が $0$ の式の配列で、そこに x を足すと、x と同じ変数を要素に持つ式の配列 y になります。

このプログラムの出力は次のとおりです:

y = {{x[0][0],0,x[0][2],x[0][3]},{0,x[1][1],x[1][2],x[1][3]},{x[2][0],x[2][1],x[2][2],x[2][3]},{x[3][0],x[3][1],x[3][2],1}}
f =
cons(x[0][0] +x[0][2] +x[0][3] == 2)
+cons(x[1][1] +x[1][2] +x[1][3] == 3)
+cons(x[2][0] +x[2][1] +x[2][2] +x[2][3] == 2)
+cons(1 +x[3][0] +x[3][1] +x[3][2] == 2)
+cons(x[0][0] +x[2][0] +x[3][0] == 2)
+cons(x[1][1] +x[2][1] +x[3][1] == 3)
+cons(x[0][2] +x[1][2] +x[2][2] +x[3][2] == 2)
+cons(1 +x[0][3] +x[1][3] +x[2][3] == 2)
x = {{1,0,1,0},{0,1,1,1},{1,1,0,0},{0,1,0,1}}

f は、方法 1 で置き換えた後の g と同じ式です。 固定した部分は、式を作る時点で定数になっています。

y には固定した値が定数として入っているので、sol(y) で行列全体をそのまま読めます。 一方、固定した x[0][1] などは f に現れないので、sol(x) で読もうとするとエラーになります（解）。

2 つの方法の比較

	方法 1: replace	方法 2: 式の配列に定数を代入
固定の違う問題を解く	f はそのままで、対応表を替えて置き換えるだけ	式の配列と式を作り直す
作る式の大きさ	固定する変数も含めたすべての変数で式を作ってから小さくする	初めから固定後の小さな式を作る
解の読み出し	qbpp::Sol(f).set(sol, ml) で全体の解に戻して読む	sol(y) で固定値も含めてそのまま読める

固定する値が違う問題をいくつも解くなら、式を 1 回作るだけで済む方法 1 が便利です。

固定する変数が多いときは、方法 2 のほうが式を作る時間もメモリも少なくて済みます。 方法 1 は、固定する変数も含めた大きな式をいったん作るためです。

なお、固定の結果として式から打ち消し合って消えた変数があると、方法 2 の sol(y) はその要素を読めずエラーになります。 方法 1 の qbpp::Sol(f).set(sol, ml) では、そうした変数は $0$ として入ります。

否定リテラルを使う式

式に否定リテラル ~x[i] を使う場合は、どちらの方法にも注意が要ります:

• 方法 1: qbpp::replace は x[i] と ~x[i] を別のキーとして扱います。 置き換える前に simplify_as_binary() を呼ぶと、2 次以下の項の ~x[i] は 1 - x[i] に展開されるので、x[i] だけをキーにすれば済みます。 3 次以上の項の否定リテラルは展開されずに残るので、その場合は ~x[i] もキーに入れます。
• 方法 2: ~ は変数にしか付けられないので、否定の式の配列は ~x から別に作り、同じ位置に反対の値（$1$ なら $0$）を代入します。

どちらも書き忘れるとエラーにはならず、固定したはずの変数が残った式になります。

番兵: 配列の外側に定数を置く

隣り合う要素の組を見る式では、配列の端の要素には片側にしか相手がありません。 配列の外側に固定値（番兵）を置くと、端の要素も内側と同じ式で扱えます。 番兵も固定した値なので、同じ 2 つの方法で置けます。

例として、最大部分配列問題を解きます。 整数の列 $w_0, \ldots, w_{n-1}$ から連続する区間を 1 つ選び、その和を最大にする問題です。 $x_i = 1$ のとき $w_i$ を選ぶとして、両端に番兵 $0$ を置いた列

\[(y_0, y_1, \ldots, y_n, y_{n+1}) = (0, x_0, \ldots, x_{n-1}, 0)\]

に対して

\[t = \sum_{i=0}^{n} (y_i - y_{i+1})^2\]

とします。各項は隣り合う 2 つの値が違うときだけ $1$ になるので、$t$ は 0 と 1 の境目の数です。 両端が $0$ なので、1 の塊が 1 つなら境目はちょうど 2 つです。 したがって、$t = 2$ の制約のもとで $\sum_i w_i x_i$ を最大にすればよいことになります。 番兵がなければ、端にある塊の境目を数えるために、$x_0$ と $x_{n-1}$ を別に足す必要があります。

$w = (-2, 1, -3, 4, -1, 2, 1, -5, 4)$ なら、答えは区間 $(4, -1, 2, 1)$ で、和は $6$ です。

方法 1: 番兵の分も変数を用意して replace で固定する

方法 1 では、番兵の分も含めて長さ $n+2$ の変数の配列 z を用意し、両端の z[0] と z[n + 1] を $0$ に固定します:

#include <qbpp/qbpp.hpp>
#include <qbpp/exhaustive_solver.hpp>

int main() {
  auto w = qbpp::array({-2, 1, -3, 4, -1, 2, 1, -5, 4});
  const size_t n = w.size();
  auto z = qbpp::var("z", n + 2);
  auto x = z(qbpp::slice(1, n + 1));

  qbpp::Expr t = 0;
  for (size_t i = 0; i <= n; ++i) t += qbpp::sqr(z[i] - z[i + 1]);
  auto f = -qbpp::sum(w * x) + 10 * qbpp::cons(t == 2);

  qbpp::MapList ml({{z[0], 0}, {z[n + 1], 0}});
  auto g = qbpp::replace(f, ml);
  g.simplify_as_binary();

  auto sol = qbpp::ExhaustiveSolver(g).search();
  auto full_sol = qbpp::Sol(f).set(sol, ml);
  std::cout << "sum = " << -sol.energy() << std::endl;
  std::cout << "x = " << full_sol(x) << std::endl;
}

このプログラムの出力は次のとおりです:

sum = 6
x = {0,0,0,1,1,1,1,0,0}

$x_i$ にあたる変数は z[i + 1] で、添字が 1 つずれます。 長さ $n$ の w と掛けるには、内側の x = z(qbpp::slice(1, n + 1)) を取り出す必要があります。 番兵の z[0] と z[n + 1] は解に含まれないので、qbpp::Sol(f).set(sol, ml) で全体の解に戻してから、内側の x を読みます。

方法 2: 番兵を含む式の配列を作る

方法 2 では、変数の配列 x の両端に定数 0 を連結した式の配列 y を作ります:

#include <qbpp/qbpp.hpp>
#include <qbpp/exhaustive_solver.hpp>

int main() {
  auto w = qbpp::array({-2, 1, -3, 4, -1, 2, 1, -5, 4});
  const size_t n = w.size();
  auto x = qbpp::var("x", n);
  auto y = qbpp::concat(0, qbpp::concat(x, 0));
  std::cout << "y = " << y << std::endl;

  qbpp::Expr t = 0;
  for (size_t i = 0; i <= n; ++i) t += qbpp::sqr(y[i] - y[i + 1]);
  auto f = -qbpp::sum(w * x) + 10 * qbpp::cons(t == 2);
  f.simplify_as_binary();

  auto sol = qbpp::ExhaustiveSolver(f).search();
  std::cout << "sum = " << -sol.energy() << std::endl;
  std::cout << "x = " << sol(x) << std::endl;
}

このプログラムの出力は次のとおりです:

y = {0,x[0],x[1],x[2],x[3],x[4],x[5],x[6],x[7],x[8],0}
sum = 6
x = {0,0,0,1,1,1,1,0,0}

変数の配列 x の添字は $x_i$ の添字のままです。 w * x はそのまま書け、解も sol(x) で読めます。 番兵は式の配列 y にだけあり、式を作るときには定数として扱われます。

番兵の 2 つの方法の比較

	方法 1: replace	方法 2: 式の配列
変数の添字	番兵の分だけずれる（$x_i$ は z[i + 1]）	$x_i$ の添字のまま
データとの演算	内側を取り出してから掛ける	w * x とそのまま書ける
解の読み出し	qbpp::Sol(f).set(sol, ml) で戻し、内側を読む	sol(x) でそのまま読める

番兵は配列の端にあるだけなので、式の大きさはどちらでもほとんど変わりません。 番兵を置く場合は、添字がずれず、データや解をそのまま扱える方法 2 が書きやすくなります。

ハーフトーン化のプログラムでは、2 次元の画像の周囲に方法 2 でゼロの番兵を置いています。